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PREFACE 


This  report  is  an  edited  version  of  the  dissertation  submitted  by  Terry  L.  Fore¬ 
man  in  partial  fulfillment  of  the  requirements  for  a  Ph.D.  degree  in  physics  at  The 
University  of  Texas  at  Austin.  Portions  of  this  work  have  been  submitted  for  pub¬ 
lication  in  the  Journal  of  the  Acoustical  Society  of  America.  This  effort  was  sup¬ 
ported  by  the  Office  of  Naval  Research  through  Contracts  N00014-80-C-0490  and 
N00014-87-K-0346. 


I.  INTRODUCTION 


Geometrical  ray  theory  propagation  models  have  consistently  been  among  the 
most  versatile  and  frequently  used  instruments  of  analysis  in  the  underwater  acous¬ 
tician’s  inventory.  The  popularity  of  these  models  derives  in  part  from  certain  com¬ 
putational  advantages,  particularly  the  ability  to  treat  range  dependent  ocean  envi¬ 
ronments  and  broadband  sources  with  relative  ease.  But  the  success  of  the  models  is 
explained  at  least  as  much  by  their  appeal  to  physical  intuition — the  easily  concep¬ 
tualized  representation  of  the  acoustic  field  by  interacting  rays  of  energy.  Decades 
of  experience  with  ray  theory  have  brought  a  sense  of  familiarity  with  rays  and  how 
they  behave.  The  behavior  of  the  acoustic  field  in  a  complicated  ocean  environment 
can  often  be  understood  by  studying  a  ray  diagram,  which  may  reveal  the  locations 
of  convergence  zones  and  shadow  zones,  regions  where  bottom  interaction  is  impor¬ 
tant  or  unimportant,  ducts  where  sound  can  propagate  long  distances  with  little 
attenuation,  and  myriad  other  details  of  the  acoustical  scene.  Despite  its  well  known 
shortcomings  at  treating  caustics,  bottom  interaction,  and  weakly  ducted  propaga¬ 
tion,  ray  theory  provides  much  of  the  conceptual  framework  on  which  ideas  about 
acoustic  propagation  are  based.  Physical  intuition  about  the  nature  and  behavior  of 
rays  has  become  quite  ingrained. 

Ray  theory  is  at  its  best  when  used  to  track  propagating  discontinuities  from 
impulsive  sources.  It  is  unsurpassed  at  o>mputing  arrival  times  and  angles  and  am¬ 
plitudes  of  disturbances  as  they  advance  into  a  quiescent  medium.  It  is  also  quite 
good  at  high  frequency  acoustic  field  predictions  for  time  harmonic  sources,  particu¬ 
larly  when  propagation  distances  are  short.  These  have  been  the  traditional  domains 
of  ray  theory  application.  In  recent  years  ray  theory  hM  been  used  increasingly  to 
model  long  range,  low  frequency  propagation,  sometimes  even  in  shallow  water,  and 
its  weaknesses  have  begun  to  show.  Yet  the  computational  power  and  convenience 
of  ray  models  and  the  insight  they  provide  encourage  their  use  even  under  these  less 
than  favorable  circumstances. 

The  Helmholtz  equation  accurately  describes  the  acoustic  field  in  the  farfield  of 
a  low  frequency  time  harmonic  source.  The  application  of  classical  ray  theory  to  the 
Helmholtz  equation  is  based  on  a  well  known  approximation.  This  dissertation  is  a 
study  of  an  idealized  ray  theory  in  which  the  approximation  has  been  circumvented. 
When  the  acoustic  field  can  be  computed  accurately  by  some  independent  means,  it 
then  becomes  possible  to  trace  out  the  rays  which  would  have  resulted  had  the  ray 
theory  approximation  not  been  made.  The  resulting  exact  ray  diagrams  provide  new 
insight  into  the  behavior  of  the  acoustic  field,  and  into  the  nature  of  classical  ray 
theory  and  the  ray  theory  approximation.  This  new  ray  theory  therefore  serves,  not 
as  a  computational  method,  but  tis  a  new  method  for  representing  and  displaying  the 
acoustic  field. 

But  the  surprising  outcome  of  numerical  experiments  is  that  exact  ray  diagrams 
bear  little  resemblance  to  their  classical  counterparts,  even  at  high  frequencies.  Some 
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of  the  more  significant  differences  are  the  following.  (1)  The  ray  path  trajectories 
depend  on  the  source  frequency  and  configuration,  and  on  the  boundaries  or  radiation 
conditions.  This  is  consistent  with  a  geometrical  representation  of  the  solution  to  the 
Helmholtz  equation,  which  is  an  elliptical  partial  differential  equation.  (2)  The  exact 
rays  intrude  into  shadow  zones  impenetrable  by  classical  rays.  The  infiltration  of 
shadow  zones  is  accomplished  naturally,  as  a  consequence  of  the  faet  that  the  field  in 
the  shadow  zone,  though  small,  is  not  nonexistent;  the  rays  do  not  take  on  complex 
parameters  in  order  to  reach  the  shadow  zone,  «is  in  some  modified  classical  ray 
theories.  (3)  Exact  ray  theory  predicts  finite  fields  at  caustics.  (4)  Last,  and  most 
significant,  the  exact  rays  never  exhibit  the  multipathing  which  is  the  hallmark  of 
classical  ray  theory.  This  item  is  closely  related  to  the  third.  If  it  were  possible 
for  the  exact  rays  to  cross  each  other  at  all  then  there  would  inevitably  be  formed 
limiting  surfaces — envelopes  to  the  rays — where  the  rays  cross  tangentially.  On  these 
surfaces  the  ray  tubes  would  pinch  off  and  the  transport  equation  would  then  predict 
nonphysical  singularities  in  the  acoustic  field. 

The  development  of  exact  ray  theory  begins  in  the  next  chapter  with  a  brief 
derivation  of  the  linear  equations  of  oceanic  acoustic  wave  propagation,  primarily  to 
establish  the  notation  and  conventions  used  in  underwater  acoustics  and  to  identify 
the  assumptions,  simplifications,  and  approximations  which  underlie  the  work  to 
follow.  Since  none  of  this  material  is  really  new,  it  seems  appropriate  to  expedite  the 
development  by  employing  the  somewhat  abstract  Lagrangian  formalism  for  classical 
fields,  as  recommended  by  Morse  and  Feshbach,[l]  to  efficiently  recapitulate  and 
consolidate  results  obtained  earlier  by  a  more  lengthy  physical  analysis.  (Boyles’ 
first  chapter  contains  a  comprehensive  development. [2]) 

Chapter  III  contains  a  development  and  reassessment  of  classical  ray  theory  which 
leads  to  the  formulation  of  an  exact  ray  theory.  Hints  begin  to  appear  that  the  ray 
theory  approximation,  far  from  being  an  innocuous  approximation  of  a  small  pertur¬ 
bative  nature,  is  actually  a  quite  drastic  assumption.  To  invoke  this  approximation 
is  to  completely  change  the  character  of  the  ray  paths. 

A  computational  method  is  presented  in  Chapter  IV  which  enables  one  to  trace 
rays  without  resort  to  the  ray  theory  approximation,  provided  a  solution  to  the 
Helmholtz  equation  is  already  available.  In  other  words,  given  a  solution  to  the 
Helmholtz  equation,  the  exact  rays  for  that  case  can  be  computed  parasitically.  The 
remainder  of  Chapter  IV  is  devoted  to  the  exploration  of  the  behavior  of  exact  rays. 
Ray  diagrams  are  constructed  for  several  cases  and  contrasted  with  the  cl2Lssical  ray 
diagrams.  Here  the  consequences  of  the  ray  approximation  become  apparent.  The 
differences  between  classical  and  exact  ray  theory  are  quite  remarkable  and  much  of 
the  chapter  is  devoted  to  explaining  the  differences. 


II.  OCEANIC  ACOUSTIC  WAVE  PROPAGATION 


We  seek  equations  which  describe  acoustic  propagation  in  a  compressible,  station¬ 
ary  fluid.  This  is  to  be  accomplished  by  determining  a  Lagrangian  density  for  the 
fluid  such  that,  when  the  time  integral  of  the  total  Lagrangian  is  minimized  according 
to  Hamilton’s  principle,  the  equations  of  motion  result.  After  imposing  simplifying 
assumptions  and  approximations  appropriate  to  the  requirements  of  this  dissertation, 
these  equations  will  turn  out  to  take  the  form  of  linear,  scalar  wave  equations  and 
Helmholtz  equations. 

A.  THE  LAGRANGIAN  DENSITY 


The  Lagrangian  density  £  of  a  conservative  system  is  the  kinetic  potential  T  -  V, 
where  T  is  the  kinetic  energy  density  and  V  the  potential  energy  density.  In  the  small 
signal  approximation,  the  kinetic  energy  density  of  the  fluid  is  simply 

^  =  |/>oV-v,  (2-1) 

where  po  is  the  ambient  density  of  the  undisturbed  fluid  and  v  is  the  particle"  velocity 
of  an  element  of  fluid,  or  fluid  particle. 

Finding  the  potential  energy  density  is  slightly  more  complicated.  In  e>rder  to 
proceed  we  use  the  fact  that  acoustic  disturbances  are  perturbations  of  the"  medium 
from  its  undisturbed  state.  Accordingly,  we  express  the  total  pressure  p'  as  the  sum 
of  the  ambient  pressure  po  and  the  excess  pressure  p  due  to  the  acoustic  perturbation: 

p'(r,0  =  po(r)  -l-p(r,<); 

but  (in  a  somewhat  inconsistent  notation)  we  write  the  total  density  p  in  terms  of 
the  ambient  density  and  the  fractional  change  in  density,  a,  due  to  the  acoustic 
disturbance: 

p{r,t)  =  po(r)[l  +  <7(r,  01-  (2-2) 

(We  are  now  following  Morse  and  Feshbach,(l]  pp.  307-309.)  When  we  spi'ak  of  an 
acoustic  disturbance  it  is  understood  that  cr  <C  1  and  that  the  acoustic  pn'ssure  is 
much  smaller  than  the  bulk  modulus  (p  <C  poc^)-  Notice  that  the  ambient  i^ressure 
and  density,  as  written,  may  depend  on  position  r  but  not  on  time  t.  We  are  as¬ 
suming  a  stationary  medium  for  the  purposes  of  this  dissertation  without  claiming 
that  temporal  fluctuations  are  unimportant.  In  fact,  the  effects  of  fluctuations  of  the 
medium  on  acoustic  propagation  are  a  very  active  area  of  research. [3 ,4] 

The  potential  energy  stored  in  an  infinitesimal  volume  of  fluid  is  equ;d  to  the 
work  done  to  compress  the  fluid  from  its  original  volume  Vo  to  its  actual  volume 
Vq  —  AV'.  Thus  the  potential  energy  density  V  of  the  fluid  is  given  by 


As  the  fluid  is  being  compressed  its  volume  decre2ises  approximately  according  to 
V  =  Vo(l  —  <y).  With  a  change  of  integration  variable  the  expression  for  the  potential 
energy  becomes 

V  =  /  pda. 

Jo 

Evidently  we  need  the  equation  of  state  which  governs  the  dependence  of  excess 
pressure  on  density  in  order  to  proceed.  For  most  fluids,  including  seawater,  the 
excess  pressure  is  proportional  to  the  fractional  change  in  density  in  the  small  signal 
approximation: 

p  =  Poc^or, 

where  c  is  the  speed  of  sound  in  the  fluid.  Then,  with  another  change  of  integration 
variable,  we  finally  obtain 


as  an  expression  for  the  the  potential  energy  density. 

With  expressions  (2-1)  and  (2-3)  in  hand  for  the  kinetic  and  potential  energy 
densities,  we  can  write  down  a  tentative  form  for  the  Lagrangian  density: 

£  =  T- V 

But  this  Lagrangian  contains  four  field  variables:  the  pressure  p,  and  the  three 
components  of  the  particle  velocity  vector  v.  The  equations  of  motion  for  this  system 
would  be  an  unwieldy  set  of  coupled  vectc.-  equations. 

We  must  pare  down  the  number  of  field  variables.  Toward  that  end,  if  we  neglect 
vorticity  then  we  may  introduce  a  scalar  velocity  potential  tp  such  that 

v{r,t)  =  Vip,  (2-4) 

and  thereby  replace  the  vector  field  v  with  the  scalar  xp.  This  approximation  is  almost 
always  valid  in  underwater  acoustics.  Moreover,  Lighthill[5]  points  out  that,  in  linear 
theory,  v  =  is  the  irrotational  part  of  the  total  velocity  field  and  does  not  interact 
with  any  steady  rotational  flow  field.  We  will  make  the  further  assumption  that  there 
is  no  steady  flow  field.  This  approximation  is  usually  well  justified  in  underwater 
acoustics,  where  typical  current  and  eddy  flow  speeds  are  ~  1  m/s,  while  the  sound 
speed  has  a  nominal  value  of  1500  m/s. 

With  the  velocity  potential  thus  defined,  the  Lagrangian  density  becomes 


(2-5) 
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which  still  contains  two  field  variables,  namely  V’  and  P-  We  will  eliminate  p  in  favor 
of  V’l  although  some  researchers  prefer  to  work  with  the  pressure  and  hence  do  the 
opposite.  In  the  absence  of  viscosity  and  external  forces,  the  pressure  and  particle 
velocity  are  related  by  (see  Morse  and  Feshbach,{l]  pp.  151-160,  308) 


dv 

Po^  =  -Vp. 


(2-6) 


By  equating  gradients  in  Eq.  (2-6)  and  the  time  derivative  of  Eq.  (2-4),  we  see  that 


(2-7) 


to  within  a  constant,  which  is  set  to  zero.  The  Lagrangian  (2-5)  can  nov;  be  written 
in  the  desired  form 


(2-8) 

with  V*  as  the  only  remaining  field  variable. 

B.  THE  ACOUSTIC  WAVE  EQUATION 

According  to  Hamilton’s  principle  for  a  conservative  system,  the  system  will 
move  so  that  the  time  average  of  the  kinetic  potential  is  an  extremum  (usually 
a  minimum).  The  equations  of  motion  which  satisfy  Hamilton’s  principle  are  the 
Euler- Lagrange  equations.  We  are  now  prepared  to  generate  a  wave  equation  by 
substituting  the  Lagrangian  density  (2-8)  in  the  Euler- Lagrange  equation 


(2-8) 


-( 

dxi  \ 


n 


The  Lagrangian  is  to  be  regarded  as  an  explicit  function  of  time,  the  coordinates,  the 
field  variable,  and  the  derivatives  of  the  field  variable  with  respect  to  time  and  the 
coordinates:  C  =  £{t,Xi,il;,dtl’/dt,drpldxi).  The  summation  convention  on  repeated 
indices  is  assumed,  so  that,  for  example,  the  particle  velocity  can  be  expressed  as 

V  =  VV’ 


f  d  d  d\  , 


and  can  be  written  as  d^t()fdxidxi.  In  this  notation  the  Lagrangian  (2-8)  takes 
the  form 

r  =  62. 

2  dxi  dxi  \dt  j 
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Then  the  equation  of  motion  turns  out  to  be 


_  £(  dC  \  _d_  /  dC  \ 

“  U(s)i 

_  j9  ^ 

5t  \  dt  j  dxi  \^dxi) 

_  /  5^0  1  dpo  50  1  5^0'\ 

\5xi5i,  po  5x,  5x,  dt'^  )  ’ 


or 


vv  +  ^.V0  =  ^ 


12-9) 


Po  c*  5<^ 

which  is  the  desired  wave  equation.  Because  of  the  presence  of  the  density  gradient 
term,  Eq.  (2-9)  will  be  referred  to  as  the  modified  wave  equation.' 


C.  THE  HELMHOL1  Z  EQUATION 


The  Fourier  transform  ^(r,u>)  of  0(r,  <)  is  defined  by 


Because  the  modified  wave  equation  (2-9)  is  linear,  0  will  satisfy  the  Fourier  trans¬ 
form  of  Eq.  (2-9),  namely 

+  k^q,  =  0,  (2-10) 

Po 

which  is  the  modified  Helmholtz  equation.  The  wave  number  k  is  defined  by 

u; 

k  =  -. 
c 

Equation  (2-10)  can  be  rendered  in  a  more  familiar  form  which  lacks  the  density 
gradient  term  by  defining  a  new  field  variable: 

(2-11) 

After  substituting  the  expression  (2-11)  for  q  in  Eq.  (2-10)  we  get  the  Helmholtz 
equation^ 

V'^q  +  K^q  =  0,  (2-12) 

'One  might  object  ihat  the  derivation  of  Eq.  (2-9)  is  inconsistent  in  that  external  forces,  here 
neglected,  must  be  present  in  order  to  create  a  density  gradient.  Indeed,  if  the  external  force  F  is  ob¬ 
tained  from  a  potential  U  by  F  =  —  VK,  then  Eq.  (2-7)  should  be  modified  to  p  =  —{podij>/dt-\-V). 
But  it  is  easily  verified  that,  for  time  independent  forces,  the  resulting  extra  terms  in  the  Lagrangian 
density  have  no  effect  on  the  wave  equation. 

^Universal  agreement  on  equation  names  is  lacking,  as  usual.  The  Helmholtz  equation  (2-12) 
is  also  called  the  reduced  wave  equation.  Morse  and  Feshbach,  whose  development  we  had  been 
following,  would  call  Eq.  (2-12)  a  Helmholtz  equation  if  k  were  a  constant,  and  a  (time  independent) 
Schrodinger  equation  otherwise.  The  distinction  is  crucial  to  their  cataloging  of  the  coordinate 
systems  in  which  the  Helmholtz  equation  (in  their  sense)  is  separable.  That  distinction  is  irrelevant 
here  and  we  refer  to  Eq.  (2-12)  as  a  Helmholtz  equation  because  it  has  the  right  form. 
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whprf> 

2  ,2,lVpoVpo  IV^po 

K  =  k  +  - - - - . 

4  po  Pq  ^  Po 

However,  it  is  rarely  necessary  in  practical  applications  to  make  this  transfor¬ 
mation;  the  density  gradient  term  in  Eq.  (2-10)  can  usually  be  ignored.  Rutherford 
and  Hawker[6]  investigated  the  sensitivity  to  variations  in  the  density  gradient  of 
solutions  to  the  depth  separated  modified  Helmholtz  equation 


(Pu 


1  dpo  du 
Po  dz  dz 


-h  {k^  -  7*)u  =  0, 


(2-13) 


where  7^  is  the  separation  constant.  They  conducted  the  sensitivity  study  for  various 
sediment  layers  underlying  the  ocean  bottom  and  found  that  solutions  to  Eq.  (2-13) 
were  quite  insensitive  to  the  details  of  the  density  profile  within  sediment  layers;  they 
were,  however,  sensitive  to  the  ratios  of  the  densities  of  two  sediment  layers  at  their 
common  interface.  Thus  they  determined  that  the  density  gradient  term  could  be 
dropped  from  Eq.  (2-13)  provided  that  the  correct  density  ratios  at  interfaces  were 
retained.  Moreover,  typical  density  gradients  in  the  sediments  often  exceed  those 
in  the  water  by  an  order  of  magnitude;  if  the  density  gradient  term  is  ignorable  in 
sediments  then  a  fortiori  it  is  ignorable  in  the  water  column. 


1.  Sources  and  boundary  conditions 

We  will  only  be  concerned  with  point,  time  harmonic  sources  with  time  de¬ 
pendence  Usually,  we  will  be  concerned  with  finding  the  Green  function  'tc 

which  satisfies  the  inhomogeneous  modified  Helmholtz  equation 

VHg  +  k'^^G  =  -<!>(r  -  To) 


for  a  source  located  at  Tq. 

At  propagation  ranges  many  times  the  depth  of  the  ocean,  acoustic  energy  is 
effectively  trapped  in  the  channel  formed  by  the  surface  and  bottom  of  the  ocean,  and 
the  ocean  begins  to  take  on  the  properties  of  a  waveguide.  As  is  often  the  case  with 
propagation  in  a  waveguide,  rough  periodicities  in  the  acoustic  field  are  frequently 
evident.  Cycle  distances  of  50  km  or  so  are  typical.  Sound  waves  are  refracted  by 
inhomogeneities  in  the  medium,  and  the  energy  propagation  paths  predicted  by  ray 
theory  become  very  complicated  at  ranges  beyond  a  cycle  distance.  Failures  of  the 
ray  theory  approximation  begin  to  manifest  at  ranges  beyond  a  cycle  distance  and 
at  frequencies  low  enough  (below  5000  Hz)  to  propagate  to  such  distances  without 
being  completely  absorbed. 

The  focus  of  this  investigation  is  on  the  exact  ray  theory  developed  in  Chap¬ 
ters  III  and  IV  rather  than  the  problems  of  modeling  acoustic  propagation  in  difficult 
ocean  environments.  It  is  desirable  to  illustrate  the  theory  using  the  simplest  possible 
idealizations  of  the  oceanic  waveguide.  With  this  in  mind,  the  simplified  boundary 
and  radiation  conditions  described  below  will  usually  suffice. 
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The  atmosphere  above  the  sea  surface  is  so  tenuous  compared  to  the  water 
that  very  little  acoustic  energy  passes  from  the  water  into  the  air.  It  is  sufficient  to 
idealize  the  air-water  interface  as  a  planar  pressure  release  perfect  reflector  at  which 
the  acoustic  field  vanishes.  One  should  note,  however,  that  considerable  research  is 
underway  on  the  effects  of  scattering  from  surface  waves  and  polar  ice  caps. 

Similarly,  even  though  acoustic  interaction  with  the  ocean  bottom  is  also 
an  important  area  of  research,  it  suffices  here  to  treat  the  bottom  as  a  plane  lying 
parallel  to  the  surface. 

The  acoustic  field  is  considered  to  consist  entirely  of  outgoing  waves.  In  a 
horizontally  stratified  waveguide  this  leads  to  the  requirement  that  the  acoustic  field 
satisfy  the  Sommerfeld  rafiiation  condition 


in  cylindrical  coordinates  where  r  is  horizontal  range  along  the  waveguide  and  z  is 
depth  below  the  surface. 

2.  Attenuation 

It  is  convenient  to  define  a  complex  velocity  potential  for  a  time  harmonic 
source  by 

V?  =  (2-14) 

where 

V*  =  Rey>.  (2-15) 

It  can  be  shown  by  substitution  that  satisfies  the  acoustic  wave  equation  (2-9)  for 

sources  with  time  dependence  Moreover,  ^  also  satisfies  the  modified  Helmholtz 
equation: 

VV  +  — -Vcp  4-  k\  =  0;  (2-16) 

Po 

in  fact,  the  modified  Helmholtz  equation  can  be  derived  by  substituting  ip  in  the 
modified  wave  equation. 

We  now  exploit  the  property  of  the  Lagrange  formulation  of  classical  fields 
that  all  Lagrangian  densities  which  generate  the  same  equations  of  motion  are  equiv¬ 
alent^  (though  not  all  Lagrangians  have  the  same  physical  content).  Here,  a  La¬ 
grangian  density  which  generates  the  modified  Helmholtz  equation  (2-16)  is 

C  =  pQ{V<f’Vp‘ -  (2-17) 

^Actually,  this  statement  is  strictly  true  only  when  the  Lagrangian  is  invariant  to  coordinate 
transformations.  We  avoid  these  complications  by  staying  with  Cartesian  coordinates  while  working 
out  the  Ruler-Lagrange  equations.  Once  the  equations  of  motion  have  been  derived  they  can  be 
rewritten  for  any  desired  coordinate  system. 


where  there  are  now  two  independent  field  variables,  tp  and  the  conjugate  field  p*. 
Substitution  of  the  Lagrangian  (2-17)  in  the  Euler- Lagrange  equation,  with  p*  as 
the  field  variable,  produces  the  modified  Helmholtz  equation 

VV  -f  — +  k^p  =  0. 

Po 

Similarly,  the  Euler- Lagrange  equation  with  respect  to  p  is 

+  Yfl.Vp^  +  itV*  =  0, 

Po 

and  the  conjugate  field  p*  is  thus  seen  also  to  satisfy  the  modified  Helmholtz  equation. 

Only  for  conservative  systems  will  the  Lagrange  function  density  take  the 
form  T  —  V.  Having  thus  excluded  dissipative  mechanisms  at  the  outset  of  our 
development  by  using  the  Lagrangian  formulation,  we  will  now  readmit  them  by  a 
subterfuge.  The  trick  is  to  add  an  extra  term  to  the  Lagrangian  (2-17): 

C  =  po  Vp-Vp* -k^pp" +  ^a^p*^-p^'^  , 

which  will  cause  the  conjugate  field  p*  to  gain  in  amplitude  as  rapidly  as  the  field  p 
attenuates,  thereby  preserving  the  total  energy  of  the  system  as  a  constant  and  thus 
permitting  the  continued  use  of  the  Lagrangian  formulation  (it  is  understood  that 
these  gains  and  losses  take  place  in  space,  not  time).  The  Euler- Lagrange  equations 
are 

VV  +  — -Vv?  +  itV  -  =  0  (2-18) 

Po  ot 

and 

VV*  +  — -h  ibV  +  =  0.  (2-19) 

Po 

Since  p  =  the  field  equations  (2-18)  and  (2-19)  become 

-f  -f-  (it^  +  iau;)<lf  =  0  (2-20) 

po 

and 

-I-  {e  -  mu;)'!'*  =  0,  (2-21) 

Po 

now  in  terms  of  and  ^'*  again.  Only  the  solution  to  Eq.  (2-20)  is  of  direct  physical 
interest,  for will  suffer  the  desired  attenuation  while  'I'*,  the  solution  to  Eq.  (2-21), 
will  grow  to  compensate. 

Morse  and  Feshbach  recommend  this  ad  hoc  device  for  introducing  dissipation 
into  the  Lagrangian  only  when  the  loss  mechanisms  are  not  well  understood.  The 
mechanisr?7s  of  altemialioii  in  the  ocean  are  actually  understood  reasonably  well[7) 
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but  are  aJso  rather  complicated;  to  treat  them  thoroughly  would  take  us  too  far 
afield.  Instead  we  embrace  equation  (2-20),  which  is  commonly  modified  to 


+  {k  +  =  0 

Po 


(2-22) 


(valid  for  a  with  a  —  jOc),  as  a  successful  phenomenological  model  of  time 
harmonic  propagation  in  a  slightly  lossy  ocean.  The  attenuation  a(r,u;)  is  usually 
determined  empirically  and  incorporates  the  effects  of  thermoviscous  and  relaxation 
loss  mechanisms  (absorption)  as  well  as  bulk  scattering  losses.  It  has  dimensions  of 
nepers  per  unit  length. 

The  remainder  of  this  dissertation  is  chiefly  concerned  with  constructing  exact 
ray  diagrams  for  various  solutions  to  the  lossy  modified  Helmholtz  equation  (2-22) 
and  with  using  exact  ray  theory  to  exhibit  the  behavior  of  the  acoustic  field. 


III.  REVIEW  OF  RAY  THEORY 


This  chapter  begins  with  what  is  essentially  the  Sommerfeld  development  of 
ray  theory, [8]  perhaps  the  simplest  of  several  alternative  derivations.  Moreover,  the 
Sommerfeld  development  directly  addresses  the  application  of  ray  theory  to  the  Helm¬ 
holtz  equation  and  it  requires  only  one  approximation,  which  is  clearly  identified.  The 
methodology  and  order  of  presentation  given  here  are  somewhat  unconventional  in 
order  to  advance  a  view  of  classical  ray  acoustics  as  an  approximation  to  an  exact 
ray  theory  in  which  the  ray  trajectories  depend  on  the  configuration  and  radiating 
frequency  of  the  time  harmonic  sources.*  Ray  theory  in  the  traditional  view  is  the 
high  frequency  limit  of  an  asymptotic  expansion  in  which  the  ray  paths  arc  indepen¬ 
dent  of  frequency.  The  two  viewpoints  are  not  incompatible  and  a  proof  of  the  latter 
claim  is  also  given. 

A.  THE  RAY  PATH  AND  TRANSPORT  EQUATIONS 

We  will  obtain  the  eikonal  and  transport  equations  of  ray  acoustics  by  further 
manipulation  of  the  Lagrangian,  and  from  these  we  will  derive  the  ray  path  equation 
and  an  energy  conservation  law. 

A  Lagrangian  density  (2-17)  which  yields  the  modified  Helmholtz  equation  in  a 
lossless  medium  was  found  in  Sec.  II.C.2  to  be 


which  becomes 

£  =  po(  (3-1 ) 

upon  using  the  defining  relation  (p  =  Henceforth  we  will  assume  that  the  den¬ 

sity  dependence  of  the  Helmholtz  equation  has  been  abolished,  either  by  a  change  of 
variables  or  by  simply  ignoring  it,  as  discussed  in  Sec.  II.C.  Then  the  Lagrangian  (3-1) 
may  be  replaced  by 

C  =  (3-2) 

as  one  may  verify  by  substituting  C  in  the  Euler- Lagrange  equation  and  showing 
that  the  Helmholtz  equation 

V^'4'  +  =  0  (3-3) 

is  thereby  recovered. 

A  Lagrangian  density  which  leads  to  the  ray  theoretical  formulation  of  tlie  Helm¬ 
holtz  equation  is  obtained  by  writing  ^  explicitly  in  terms  of  its  real  amplitude  A 
and  phase 

^  =  Ae'*, 

*  A  somewhat  similar  theoretical  development  by  Floyd[9]  does  not  retain  the  coiico|)(.  of  a  geo¬ 
metrical  energy  conservation  law,  which  is  a  central  feature  of  the  ray  theory  described  here. 
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and  substituting  in  the  Lagrangian  (3-2)  to  produce  a  new  Lagrangian 


£  =  {VA  +  iAV$).( VA  - 

=  VA-VA  +  A^V^‘V(^-k^A^. 

When  the  Euler- Lagrange  machinery  is  applied  to  the  new  Lagrangian  it  generates 
two  coupled  equations  for  the  new  field  variables  A  and  One  of  the  Euler- Lagrange 
equations  is 


_  d  I  dC  \  dC 

a  A 

=  2( 


aM 

dxidxi 


-  -I- ArM 


)■ 


or 


A 


(3-4) 


Equation  (3-4)  will  be  referred  to  as  the  exact  eikonal  equation  in  order  to  distinguish 
it  from  the  more  familiar  eikonal  equation  of  classical  ray  acoustics,  which  will  appear 
later.  Likewise,  $  will  be  called  the  exact  eikonal. 

The  other  Euler- Lagrange  equation  is 


^  d  /  a£  \ 

"  V^/ 

=  2(: 


d± 

dx, 

dA^ 

dxi  dxi 


dC 


dxidx 


-)■ 


or 


V4 

2__.V$  -f  =  0, 
A 


(3-5) 


which  is  the  transport  equation. 

The  eikonal  equation  leads  to  an  equation  for  the  paths  of  acoustic  energy  prop¬ 
agation,  or  ray  paths.  The  transport  equation  lejwls  to  a  geometrical  energy  conser¬ 
vation  law  in  a  form  useful  for  the  calculation  of  the  sound  intensity  at  any  point 
along  a  ray  path.  The  eikonal  and  transport  equations  are  coupled,  first  order  in 
Vy4/A  and  V$,  and  nonlinear;  they  are  entirely  equivalent  to  the  second  order,  lin¬ 
ear  Helmholtz  equation  (3-3).  They  may  also  be  obtained,  <is  in  the  Sommerfeld 
development,  by  direct  substitution  of  i4e’*  in  the  Helmholtz  equation. 

To  bring  the  concept  of  rays  into  the  development,  first  define  the  unit  vector  t 
to  be  the  outward  normal  to  a  surface  of  constant  phaae  $  at  the  point  r,  and 
further  define  t  to  be  the  tangent  vector  to  the  ray  path  at  r.  (In  some  ray  theories, 
particularly  those  applicable  to  moving  media,  the  ray  paths  are  not  necessarily 
perpendicular  to  the  surfaces  of  constant  phase.)  It  is  also  convenient  to  define  the 
wave  vector  K  by 


K  =  V4>, 


so  that  K  =  Kt.  The  ray  path  trajectory  is  the  solution  to  dr/ds  =  t,  where  ds  is 
the  element  of  ray  path  length.  The  exact  eikonal  equation  then  becomes 

(3-6) 

The  transport  equation  (3-5)  may  be  expressed  as  a  vanishing  divergence: 

V*(A*K)  =  0.  (3-7) 

B.  CONSERVATION  OF  ENERGY 

Departing  slightly  from  the  usual  order  of  the  Sommerfeld  development,  we  will 
follow  up  on  the  discussion  of  the  transport  equation  before  deriving  the  ray  path 
equation.  The  formulation  of  the  transport  equation  (3-7)  as  a  vanishing  divergence 
suggests  applying  Gauss’s  theorem.  Accordingly,  let  S  be  the  closed  surface  formed 
by  a  tube  of  rays  and  the  intersection  of  the  tube  with  two  surfaces  of  constant  phase, 
5i  and  52,  as  shown  in  Fig.  3.1.  Let  V  be  the  volume  enclosed  by  5.  Then  by  Gauss’s 


Figure  3.1  Conservation  of  energy  in  ray  acoustics.  The  time  averaged  energy  flux 

THROUGH  THE  SURFACE  OF  CONSTANT  PHASE,  S\,  EQUALS  THE  FLUX  THROUGH  Sj. 


theorem, 

j^ViA‘^K)dV  =  J^A^K’udS  =  0,  (3-8) 

where  n  is  the  outward  unit  vector  normal  to  5.  The  only  contributions  to  the 
integral  are  from  the  end  caps  Si  and  S2  since  K*n  =  0  by  definition  on  llie  sides 
of  5  formed  by  the  ray  tube.  And  since  K  is  tangential  to  the  ray  path,  then  also  by 
definition. 


n  = 


— t  on  5i 
t  on  Si 


(3-9) 


Using  Eq.  (3-9)  in  Eki-  (3-8)  yields 


I 

I 

i 


or 


=  J^A^K-ndS 

=  /  A^KAdS-f  A^K-tdS, 

JSi  JSi 

f  A^KdS=  /  A^KdS, 

JSi  Js, 


(3-10) 


which  is  the  integral  form  of  the  transport  equation. 

The  significance  of  Ek^.  (3-10)  becomes  apparent  upon  identifying  -pouiA^K  as 
the  time  averaged  acoustic  intensity,  or  energy  flux.  The  intensity  of  a  classical 
field  has  dimensions  of  force  per  unit  area  times  velocity.  We  therefore  expect,  by 
dimensional  analysis,  that  the  instantaneous  intensity  J  is  given  by 


J  =  pv  (3-11) 

since  pressure  has  units  of  force  per  unit  area.  (A  complete  derivation  is  given  by 
Boyles,[2]  pp.  30-34.)  The  energy  flux  J  is  sometimes  called  the  Poynting  vector, 
although  that  term  is  more  commonly  applied  to  the  analogous  electromagnetic  field 
intensity. 

In  terms  of  the  complex  velocity  potential  the  particle  velocity  (2-4)  is 


V  =  Re 


(3-12) 


and  the  pressure  (2-7)  is 


p  =  -poRe 


(^g-iwt) 


=  — powIm(^e  *"*) 


(3-13) 


Then,  using  ^  =  Ae'*  and  Eqs.  (3-12)  and  (3-13)  for  v  and  p,  the  expression  for  the 
intensity  (3-11)  becomes 


VA 


J  =  PqU>A^  [  V$  sin*($  —  ujt) - —  sin($  —  u;t)  cos($ 


-u;<)^  . 


The  time  average  of  the  intensity,  (J),  is  obtained  by  integrating  J  over  the  period  T: 

with  T  =  u}l2ir.  The  integral  is  easily  calculated  and  the  result  takes  the  form 


(J)  =  ipou;A*K 


(3-14) 


if  one  uses  the  definition  =  K.  Thus  the  time  averaged  intensity  is  directed  along 
the  ray  path,^  although  the  instantaneous  flux  vector  points  in  the  direction  of  v, 
which  is  time  dependent  and  generally  not  parallel  to  the  ray  path. 

In  view  of  Elq.  (3-14),  the  integral  form  of  the  transport  equation  (3-10)  can  be 
written  as 


which  is 


(3-15) 


This  is  an  energy  conservation  law  in  geometrical  form.  It  states  that  all  of  the  energy 
flux  which  parses  through  the  surface  S\  also  passes  through  52.  The  existence  of 
such  a  law  helps  to  explain  the  popularity  and  versatility  of  ray  theory. 

We  have  established  not  only  that  the  ray  paths  are  physically  meaningful  in 
that  they  are  the  trajectories  of  energy  propagation  (although  we  have  not  yet  seen 
how  to  compute  the  trajectories),  but  wc  also  have  an  intuitively  appealing  picture 
of  how  the  intensity  varies  along  a  ray  path  as  the  ray  tube  expands  and  contracts. 
The  energy  conservation  law  (3-15)  also  readily  lends  itself  to  computation.  Ray 
theory  is  unique  in  its  ability  to  present  a  qualitative  and  informative  view  of  energy 
propagation  and  to  permit  straightforward  calculation  of  the  acoustic  field. 

It  should  be  emphasized  that  no  approximations  have  been  made  up  to  this 
point;  the  geometrical  energy  conservation  law  (3-15)  is  not  a  ray  theory  approxi¬ 
mation,  as  is  sometimes  stated;  it  does  not  arise  from  the  “neglect  of  diffraction,” 
or  the  “local  approximation  of  the  field  by  a  plane  wave,”  or  any  other  classical 
ray  theory  approximation.  Moreover,  Eqs.  (3-10)  and  (3-15)  are  fully  applicable  to 
finite  sources;  nothing  in  the  development  requires  that  the  ray  tubes  include  only 
infinitesimal  volumes.  We  will,  in  the  next  section,  resort  to  an  approximation  in 
order  to  be  able  to  compute  the  ray  trajectories,  but  this  approximation  is  made  in 
the  ray  path  equation,  not  the  transport  equation. 


C. 

We  have  shown  that  there  are  trajectories  along  which  sound  propagates.  An 
equation  for  these  paths  may  be  obtained  by  computing  the  gradient  of  K^: 

=  K-K; 

2KVK  =  2(K.V)K  +  2Kx(VxK). 

can  be  shown  (see  Boyle8[2])  that  the  result  (3-14)  can  also  be  found  from  (J)  =  ^  Re(PV*) 
for  time  harmonic  fields,  where  the  complex  particle  velocity  is  defined  by  V  =  =  V’i’e"*"* 

and  the  complex  pressure  is  given  by  P  =  —podipfdl  = 


THE  RAY  PATH  EQUATION 


i 


Since  V  xK  =  V  X  V$  =  0,  this  simplifies  to 


■  (r')- 

- 


±  Ik^\  = 

ds\  dsj 


(3-16) 


Equation  (3-16)  is  a  ray  path  equation;  its  solutions  are  the  energy  conserving  ray 
path  trajectories  described  in  the  previous  section. 

It  is  worth  reemphasizing  that  no  approximations  have  yet  been  made;  the  Helm¬ 
holtz  equation  has  been  recast  in  what  might  be  termed  an  exact  ray  theoretical 
formulation.  The  simultaneous  solution  of  the  ray  path  equation  (3-16)  and  the 
transport  equation  (3-7)  would  yield  A  and  $  and  thus  determine  ’P,  the  solution  to 
the  Helmholtz  equation. 

The  difficulty  is  that  solving  the  ray  path  equation  requires  knowledge  of  /f, 
which  in  turn  requires  knowledge  of  A,  as  may  be  seen  by  examining  the  eikonal 
equation  (3-6).  And,  of  course,  the  ray  paths  must  be  known  in  order  to  solve  the 
transport  equation  for  A.  We  have  come  full  circle.  In  principle,  the  exact  eikonal  and 
transport  equations  can  be  solved  simultaneously,  for  example  by  numerical  methods, 
but  that  will  usually  prove  more  difficult  than  solving  the  original  Helmholtz  equation. 

Yet  this  “exact  ray  theory”  nevertheless  proves  useful,  not  as  a  genera!  purpose 
computational  procedure,  but  as  an  investigative  and  display  technique.  When  the 
solution  to  the  Helmholtz  equation  is  available  by  independent  means  then  the  cycle 
described  in  the  last  paragraph  can  be  broken;  givtn  the  solution  to  the  Helmholtz 
equation  one  could  compute  A  and  then  K  and  proceed  to  trace  rays.  ( A  more  prac¬ 
tical  computational  method  for  tracing  exact  rays  is  explained  in  the  next  chapter.) 
Exact  ray  diagrams  are  a  new  alternative  to  propagation  loss  curves  and  inten¬ 
sity  plots  as  a  means  of  displaying  and  analyzing  the  acoustic  field.  As  we  shall 
see  in  Chapter  IV,  exact  ray  diagrams  usually  differ  strikingly  from  their  classical 
counterparts. 

The  source  of  the  difficulty  with  exact  ray  theory  as  a  computational  method 
is  that  the  eikonal  and  transport  equations  are  coupled:  one  cannot  determine  K 
exactly  without  first  knowing  A.  In  classical  ray  acoustics  the  desired  decoupling 
is  achieved  by  resorting  to  an  approximation.  If  VM/A  <  jb*,  then  by  neglecting 
V^Af A  in  the  eikonal  equation  one  obtains 

K  w  k.  (3-17) 

This  is  the  classical  ray  theory  approximation;  it  is  the  only  approximation  required 
to  obtain  classical  ray  acoustics,  but  it  has  far-reaching  consequences. 
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It  can  be  shown  that  the  ray  theory  approximation  will  fail  wherever  the  sound 
speed  gradient  changes  significantly  over  a  wavelength. [10]  In  underwater  acoustics, 
such  rapid  variations  can  occur  at  fronts  where  cold  and  warm  ocean  waters  meet, 
or  when  one  attempts  to  perform  low  frequency  ray  calculations  in  the  ocean  bottom 
sediments.  However,  the  approximation  is  also  observed  to  fail  regularly  even  when 
the  ocean  environment  harbors  no  regions  of  rapid  sound  speed  variation  at  all. 
Regions  where  the  acoustic  field  amplitude  undergoes  rapid  spatial  variations  (t.e., 
where  V^AjA  ~  and  the  approximation  is  invalid  by  definition)  are  simply  not 
always,  or  even  often,  collocated  with  or  causally  related  to  regions  of  rapid  sound 
speed  variations.  These  problematical  regions  include  the  vicinity  of  caustics  (to 
be  discussed  at  length  in  Chapter  IV)  and  ocean  surface  ducts.  Because  of  their 
ubiquity,  these  regions  are  far  more  troublesome  to  the  routine  use  of  ray  theory 
than  rare  pathological  features  in  the  sound  speed  field. 

Nevertheless,  the  ray  theory  approximation  does  decouple  the  eikonal  and  trans¬ 
port  equations  and  thus  removes  the  dependence  of  K  on  A.  The  ray  path  equa¬ 
tion  (3-16)  then  becomes 


Since  is  a  known  function  of  the  frequency  and  sound  speed  one  can  solve  the  ray 
path  equation  and  trace  rays  at  will. 

If  we  define  an  index  of  refraction  n  in  terms  of  a  constant  reference  sound 
speed^  Cq: 

n  =  cq/c, 

then  k  =  ujuIcq  and  the  ray  path  equation  becomes 


_d 

ds 


=  Vn, 


(3-18) 


which  is  perhaps  the  most  familiar  form  of  the  classical  ray  path  equation.  In  the 
ray  theory  approximation  the  transport  equation  becomes 


V- 


=  0, 


and  is  generally  restricted  to  infinitesimal  ray  tubes  for  reasons  which  will  become 
apparent  later. 

Most  fluids,  including  seawater,  are  only  very  slightly  dispersive;  that  is,  the 
dependence  of  c  or  n  on  u;  is  negligible  for  our  purposes.  The  classical  ray  paths, 
which  are  solutions  to  Eq.  (3-18),  are  consequently  independent  of  frequency.'* 

^The  natural  choice  for  the  value  of  cq  in  optics  is  the  speed  of  light  in  a  vacuum.  There  is  no 
analogous  preferred  reference  sound  speed  in  acoustics  and  many  workers  simply  use  1/c,  sometimes 
called  the  “slowness,”  in  place  of  n. 

'*When  attenuation  is  taken  into  consideration  the  wave  number  k  or,  equivalently,  the  sound 
speed  takes  on  a  small  imaginary  part  which  is  frequency  dependent.  A  lossy  medium  is  therefore 
also  a  dispersive  medium. 


U 
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By  contrast,  K  is  &  complicated  function  of  frequency  because  of  its  dependence 
on  A,  which  is  sensitive  to  frequency.  The  ray  trajectories  which  are  solutions  to 
the  exact  ray  path  equation  (3-16),  and  which  are  the  paths  of  acoustic  energy 
propagation,  are  inherently  frequency  dependent.  We  generally  expect  dispersion  in 
an  inhomogeneous  medium  even  if  the  sound  speed  itself  is  not  an  explicit  function 
of  frequency. 

D.  THE  TREATMENT  OF  ATTENUATION  IN  RAY  THEORY 

Let  us  see  how  the  inclusion  of  attenuation  modifies  the  ray  theory  analysis 
presented  earlier  for  lossless  media.  It  was  shown  in  Sec.  I1.C.2  that,  for  a  source 
with  time  dependence  c~*"‘,  one  may  account  for  the  effects  of  attenuation  to  first 
order  by  adding  a  small  imaginary  part  to  k;  that  is,  k  k  ia.  Accordingly,  we  will 
seek  a  first  order  perturbative  solution  to  the  lossy  Helmholtz  equation  by  expanding 
about  the  solution  to  the  lossless  equation  in  the  small  parameter  a,  with  c  as  an 
ordering  parameter. 

Let  ^0  =  Ae'*  denote  the  solution  to  the  lossless  Helmholtz  equation  V^'I'o  + 
jt^^o=  0  (recall  that  A  and  0  are  real).  The  lossy  Helmholtz  equation  is  written  as 

+  {k  +  =  0,  (3-19) 

and  the  solution  proposed  for  it  as 


m  (3-20) 

The  factor  e"'^  accounts  for  the  expected  spatial  damping  of  the  field  amplitude;  it 
was  written  in  this  form  in  anticipation  of  finding  that  /?  =  / ads,  which  gives  the 
usual  path  integral  correction  for  attenuation  in  ray  theory. 

After  substituting  the  expression  (3-20)  for  ^  in  Eq.  (3-19)  and  linearizing  in  t, 
one  obtains 


4-  k^% 

-  e[l3{V^^o  +  k^%)  +  2VVI>o-V^  +  -  2i%ak]  =  0. 

The  term  of  order  c®  and  the  parenthetical  expression  in  the  term  of  order  e*  simply 
reaffirm  the  Helmholtz  equation  and  vanish  identically,  leaving 

=  2i^oak.  (3-21) 

But 

__  ^  {VA  ,,dT\ 

which,  in  the  ray  theory  approximation  '^AjA  k,  simplifies  to 

V%  =  (3-22) 

as 


(3-23) 


Substitution  of  Eq.  (3-22)  for  V^o  in  Eq.  (3-21)  gives 

2iA-V0  +  =  2tak. 

as 

Since  all  of  the  variables  in  Eq.  (3-23)  are  explicitly  real,  the  real  and  imaginary 
parts  of  Eq.  (3-23)  must  vanish  separately,  leaving  a  pair  of  coupled  equations: 

W  =  0,  (3-24) 

and 

=  a.  (3-25) 

as 

By  using  the  operator  identity  dr/ds*V  =  d/ds,  Eq.  (3-25)  can  be  written  as 

d/3 


or,  upon  integrating  with  respect  to  path  length, 

rads'.  (3-26) 

Jo 

It  is  standard  practice  in  ray  calculations  to  treat  attenuation  by  applying  th('  formula 

=  Ae-^e'*, 

where  13  is  given  by  the  path  integral  (3-26).  Such  path  integrals  are  easy  to  compute 
numerically. 

Equation  (3-24)  is  essentially  an  assertion  that  /?  is  slowly  varying.  This  asser¬ 
tion  is  well  justified  in  the  ray  theory  approximation  and  in  view  of  the  path  integral 
formulation  (3-26)  for  0. 

A  similar  analysis  which  does  not  invoke  the  ray  theory  approximation  leads 
to  path  integral  corrections  for  both  the  amplitude  and  phase.  Nevertheless,  the 
corrections  still  take  the  form  of  integrals  over  the  exact  ray  paths  of  the  lossless 
Helmholtz  equation.  We  will  henceforth  treat  only  the  lossless  Helmholtz  equation 
since  corrections  for  attenuation  are  easily  computed  and  do  not  alter  the  exact  ray 
paths  (at  least  to  first  order). 

E.  REMARKS  ON  THE  NATURE  OF  CLASSICAL  RAY  ACOUSTICS 

It  was  shown  in  the  preceding  section  that  classical  ray  acoustics  can  provide  only 
approximate  solutions  to  the  Helmholtz  equation.  The  nature  of  the  approximation 
is  explored  further  in  this  section,  where  it  is  shown  that  classical  ray  acoustics  is  the 
high  frequency  limit  of  a  nonuniform  asymptotic  expansion  in  which  the  ray  paths  are 
independent  of  frequency.  Classical  ray  acoustics  can  also  provide  exact  solutions  to 
an  important  class  of  problems  related  to  propagating  discontinuities  in  the  acoustic 
field. 


20 


1.  Propagating  discontinuities  in  the  acoustic  field 

Consider  an  inhomogeneous  quiescent  medium  which  is  suddenly  disturbed 
by  activating  a  point  source,  not  necessarily  time  harmonic.  We  ask  how  much  time 
must  elapse  before  the  disturbance  can  reach  a  distant  observation  point.  Since  the 
rate  of  propagation  of  the  disturbance  at  any  point  in  the  medium  is  equal  to  the 
local  speed  of  sound  (or  else  what  does  “sound  speed”  mean  for  small  signals  in  a 
medium  at  rest?),  then  the  minimum  travel  time  must  be  given  by  the  path 
integral 


-It 

r  {dx^  +  dy^  +  dz^y!^ 

J  c{x,y,z) 

r  [{dxldsf  +  {dyidsf  +  {dzldsfyl^ 
J  c{x,y,z) 

J  c{x,y,z) 

=  fi^ds 

J  c{Xi)  ’ 


(3-27) 


where  primes  denote  differentiation  with  respect  to  path  length  along  a  trajectory 
yet  to  be  determined.  Application  of  the  variational  calculus  reveals  that  the  path 
itself  is  the  solution  to[ll,12] 


/dF\  dF 

l^axj  ax.""’ 

(3-28) 

c{xi)  • 

(3-29) 

where 


Upon  substituting  Eq.  (3-29)  for  F  in  Eq.  (3-28)  and  using  the  auxiliary  condition 
=  1,  the  path  turns  out  to  satisfy  the  classical  ray  path  equation 


In  other  words,  no  disturbance  can  reach  the  observation  point  due  to  the  activity 
of  the  source  before  a  minimum  amount  of  time  tmin  has  elapsed,  where  is  given 
by  the  classical  ray  path  integral  (3-27).* 

® Nevertheless,  propagation  times  apparently  indicative  of  supersonic  transmission  speeds  do  oc¬ 
cur  occasionally  even  in  linear  theory,  particularly  during  attempts  to  construct  synthetic  lime  .seru« 
by  Fourier  synthesis  of  solutions  to  the  lossy  Helmholtz  equation.  'I'he  problem  usually  is  due  either 


A 


The  preceding  argument,  though  inspired  by  Fermat’s  principle  and  Hamil¬ 
ton’s  principle  of  least  action,  makes  i.o  reference  to  physical  processes  (in  particular, 
the  wave  equation  was  never  invoked);  hence,  it  gives  no  assurance  that  any  acous¬ 
tical  energy  will  actually  flow  along  the  classical  ray  paths  and  arrive  at  the  times 
predicited  by  ray  theory.  Nevertheless,  it  does  serve  to  establish  a  lower  bound  on 
transmission  times.  But  there  is  a  considerable  body  of  literature,  going  back  at 
least  to  ChristofFel[15]  in  1877,  which  shows  that  a  surface  across  which  the  time 
derivatives  of  the  acoustic  field  are  discontinuous  propagates  according  to  the  clas¬ 
sical  eikonal  equation.  Luneburg[16]  obtained  a  similar  result  for  electromagnetic 
propagation.  Heller[17]  obtained  a  generalized  eikonal  which  allows  for  movement  of 
the  fluid  medium.  And  Keller[l8]  showed  that  the  transport  equation  governs  the 
amplitudes  of  propagating  discontinuities.  The  conclusion  to  be  drawn  from  this  is 
that  classical  ray  theory  is  essentially  a  broadband  theory  best  suited  to  the  analysis 
of  propagating  fronts  generated  by  pulsive  sources. 

2.  Ray  theory  as  a  high  frequency  approximation 

We  seek  an  asymptotic  solution  of  the  Helmholtz  equation  for  high  frequencies 
by  expanding  in  the  small  quantity  I /a;.  Referring  to  the  Sommerfeld  representation 
of  the  acoustic  field,  ^  =:  Ae'*,  we  express  the  exact  eikonal  as 


=  u;T(r), 

while  the  amplitude  is  written  as  an  asymptotic  series: 


(3-30) 


E  (-) 


(3-31) 


Note  that  t,  which  will  determine  the  ray  trajectories,  is  independent  of  frequency. 
After  substituting  for  ^  in  the  Helmholtz  equation,  which  is  written  in  the  form 


-I- =  0, 

<r 


one  obtains 


-I-  j  ^2 VA„-Vt  -h  An ^  VM„  =  0.  (3-32) 

cal  analog  of  the  Kramers-Kronig  relations, [13, 14]  or  to  a  failure  of  a  saddle  point  approximation 
or  other  approximation  method  used  to  estimate  group  velocities.  Either  condition  can  result  in 
predictions  of  spurious  precursors  to  the  earliest  arrival  predicted  by  ray  theory.  Fortunately,  the 
intensities  predicted  for  these  premature  arrivals  are  usually  low  enough  that  they  can  be  tolerated 
without  requiring  corrective  measures.  Sometimes  these  acausal  precursors  are  erroneously  accepted 
at  face  value  and  attributed  to  mysterious  “diffracted  energy  paths.” 


Since  Eq.  (3-32)  is  an  identity  in  l/w,  the  coefficients  of  each  power  of  l/w  must 
vanish  independently.  At  lowest  order  in  1  /u>  we  recover  the  eikonal  equation 

Vt-Vt  =  *4, 

while  the  first  order  equation  turns  out  to  be  the  transport  equation 

2VAq-Vt  +  =  0. 

Thus  classical  ray  acoustics  is  established  as  the  asymptotic  solution  to  the  Helmholtz 
equation  in  the  high  frequency  limit. 

The  remaining  terms  in  Eq.  (3-32)  can  be  consolidated  in  like  powers  of  1/u; 
and  written  in  the  form 

oo  ^  ■  V  n+2  /  \ 

Ejjj  (2VA„.^,.Vr  +  W  -  =  0, 

which  establishes  the  recursion  relation 

2VA„+,.Vr  +  A„+,VV  =  (3-33) 

for  the  higher  order  terms  in  the  asymptotic  expansion.  Note,  though,  that  these 
terms  affect  only  the  amplitude  A  and  not  the  ray  paths,  which  are  independent  of 
frequency.  The  inRnitesimal  ray  tube  cross  sections  are  not  prevented  from  vanish¬ 
ing  at  caustics;  hence,  this  perturbative  treatment  continues  to  predict  nonphysical 
singular  acoustic  fields  at  caustics  to  any  order  in  the  expansion.  The  expansion  is 
thus  seen  to  be  nonuniform  when  caustics  are  present. 

Equation  (3-33)  would  seem  to  offer  a  means  of  obtaining  improved  ray  theory 
calculations  by  iteration,  at  least  in  regions  free  of  caustics.  The  ray  paths  and  the 
surfaces  of  constant  phase  form  an  orthogonal  coordinate  system.  Pitre[19]  points 
out  that  Eq.  (3-33)  can  be  rendered  as  an  inhomogeneous  transport  equation  in  ray 
coordinates,  and  then  converted  into  a  ray  path  integral  to  be  evaluated  numerically. 
Chen  and  Ludwig[20]  exploited  this  procedure  to  evaluate  the  importance  of  ray 
theory  corrections. 

However,  the  iterative  correction  scheme  has  liabilities  as  a  practical  com¬ 
putational  procedure.  The  use  of  ray  coordinates  is  not  the  problem;  ray  theory 
computer  models  which  compute  the  rate  of  geometrical  ray  path  spreading  already 
have  the  information  necessary  to  compute  the  Jacobian  of  the  transformation  from, 
say,  cylindrical  coordinates  to  ray  coordinates.  But  is  not  readily  available  and, 
in  any  case,  the  corrections  are  of  dubious  value  even  when  they  can  be  generated. 
The  corrections  are  negligible  far  from  caustics  at  the  frequencies  of  interest  in  long 
range  propagation.  But  if  caustics  are  present  then  the  corrections  permit  one  to 
approach  the  caustic  only  a  little  more  closely  before  the  ray  approximation  fails;  no 
amount  of  iteration  will  remove  the  singularity  at  the  caustic  or  substantially  improve 
the  calculations  near  it.  A  glance  ahead  at  the  propagation  loss  curves  of  Fig.  4.5  on 
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page  37  will  make  this  point  clearer.  Figure  4.5  is  a  comparison  of  the  propagation 
loss  computed  using  ray  theory  and  also  using  Ludwig’s  uniform  asymptotic  formula, 
which  is  valid  on  and  near  the  caustic.  The  ray  theory  propagation  loss  curves  are 
indistinguishable  from  the  uniform  asymptotic  propagation  loss  curves  right  up  to 
the  caustic,  where  ray  theory  becomes  singular.  The  ray  theory  corrections  offered 
by  Eq.  (3-33)  are  unnecessary  far  from  the  caustic  and  of  little  use  near  it. 


IV.  ILLUSTRATIVE  APPLICATIONS  OF  EXACT  RAY  THEORY 

Several  derivations  of  classical  ray  theory  were  presented  in  the  preceding 
chapter,  and  it  was  hinted  that  if  the  ray  theory  approximation  were  avoided  then 
the  resulting  rays  would  be  quite  different.  This  chapter  begins  with  the  development 
of  a  simple  computational  procedure  for  numerically  constructing  exact  ray  diagrams 
whenever  a  solution  to  the  Helmholtz  equation  is  available.  The  remainder  of  the 
chapter  is  given  over  to  demonstrations  of  the  method. 


Whenever  a  solution  to  the  Helmholtz  equation  is  known,  it  is  possible  to  de¬ 
termine  what  the  ray  path  trajectories  would  be  if  the  ray  theory  approximation  is 
not  invoked.  To  see  how  this  comes  about,  let  us  suppose  that  a  solution  ^  to  the 
Helmholtz  equation  is  known  for  some  particular  source  configuration  and  wave  num¬ 
ber  k{r)  in  a  lossless  medium.  After  computing  the  gradient  of  the  defining  equation 
for  A  and  namely  ^  =  Ae'*,  one  obtains 


-f  t’K. 


Since  A  and  K  =  are  real,  evidently 


¥-(")■ 


(4-1) 


The  latter  is  a  first  order  diflferential  equation  for  the  ray  path.  It  is  desirable  for 
computational  reasons  to  eliminate  K  from  the  ray  path  equation  by  defining  a  new 
path  parameter  a  such  that  K  da  =  ds.  Then  Eq.  (4-1)  becomes 


(4-2) 


One  traces  a  ray  by  integrating  Eq.  (4-2),  usually  by  numerical  methods. 

The  equivalence  of  the  first  order  ray  path  equation  (4-1)  and  the  second  order 
equation  (3-16)  can  easily  be  verified  by  substituting  Ae'*  for  'I'  in  Eq.  (4-1).  One 
recovers  the  defining  equation  for  the  ray  path  tangent  vector,  Kdrfds  =  from 
which  Eq.  (3-16)  was  derived. 

One  can  also  calculate  K  by 


I\  =  Im 


(?)!■ 


(4-3) 
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which  result  is  achieved  by  computing  K*K  using  Eiq.  (4-1).  Also  available  is  the 
junplitude  A,  by  way  of 

That  the  ray  path  equation  (4-2)  is  a  first  order  differential  equation  has  several 
remarla.ble  consequences,  two  having  to  do  with  the  computational  process  of  solving 
ray  path  equations,  or  tracing  rays,  and  another  having  to  do  with  the  qualitative 
behavior  of  the  paths  themselves.  Only  one  initial  condition  is  required  in  order  to 
solve  a  first  order  differential  equation,  while  a  second  order  equation  requires  two 
initial  or  boundary  conditions  for  its  solution.  To  trace  a  classical  ray  one  must  solve 
the  second  order  ray  equation  (3-18),  which  usually  means  that  the  ray  trace  must 
begin  at  the  source  (although  the  source  need  not  be  a  point;  it  may  be  a  surface), 
because  only  at  the  source  are  both  the  initial  position  ani  direction  of  the  ray 
known.  By  contrast,  given  any  starting  point,  the  starting  direction  of  an  exact  ray 
at  that  point  is  given  by  Eq.  (4-2)  because  it  is  first  order.  Hence,  an  exact  ray  trace 
may  commence  anywhere. 

Another  useful  consequence  of  the  exact  ray  equation  being  first  order  is  that 
eigenrays  for  arbitrary  points  can  often  be  found  by  tracing  the  ray  backwards  from 
an  observation  point  to  the  source.  This  is  accomplished  by  requiring  the  numerical 
integrator  to  solve  Eq.  (4-2)  for  decreasing,  rather  than  increasing,  values  of  a.  The 
resulting  ray  trace  will  proceed  from  the  observation  point  right  back  to  the  source. 

The  above  should  not  be  construed  to  mean  that  exact  ray  theory  is  a  computa¬ 
tional  tool  superior  to  classical  ray  theory.  On  the  contrary,  classical  ray  theory  can 
generate  approximate  solutions  to  the  Helmholtz  equation  ab  initio,  whereas  exact 
ray  tracing  requires  that  a  solution  already  be  known.  Nevertheless,  the  facility  with 
which  exact  rays  can  be  traced  forwards  or  backwards  from  any  starting  position 
makes  exact  ray  theory  a  versatile  and  illuminating  alternative  means  of  displaying 
and  analyzing  the  acoustic  field. 

As  promised  earlier,  we  can  also  glean  an  insight  into  the  qualitative  behavior  of 
exact  rays  from  the  fact  that  the  ray  path  equation  is  first  order.  Since  the  particle 
velocity  is  a  me<isurable  quantity  and  must  therefore  be  single  valued  at  any  position 
r,  and  since  the  velocity  potential  satisfies  the  acoustic  Helmholtz  equation,  it  follows 
that  the  velocity  potential  is  single  valued  and  twice  differentiable.  Hence,  the  right 
hand  side  of  Eq.  (4-2)  is  single  valued  at  r.  But  this  self-evident  and  seemingly 
innocuous  observation  means  that  the  vector  dr /da  tangent  to  the  ray  path  is  unique, 
which  in  turn  forces  us  to  conclude  that  only  one  ray  passes  through  each  point  r. 
There  is  exactly  one  eigenray  from  a  source  to  an  observation  point.  The  uniqueness 
of  exact  eigenrays  stands  in  stark  contrast  to  classical  rays,  which  characteristically 
exhibit  multipathing.  This  is  an  important  realization  and  a  source  of  considerable 
confusion;  we  will  return  to  it  several  times. 

The  remainder  of  this  chapter  is  devoted  to  demonstrations  of  exact  ray  analysis. 
Each  case  selected  for  analysis  was  chosen  quite  frankly  for  its  pedagogical  value  by 
the  following  criteria;  (1)  Each  case  should  serve  as  a  canonical  example  of  some 
problem  frequently  encountered  in  acoustics.  For  example,  the  Kormilitsin  problem 
has  as  its  dominant  feature  a  smooth,  convex  caustic.  The  “harmonic  oscillator” 
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waveguide  was  selected  as  a  model  of  propagation  in  a  duct.  The  circular  piston  is 
famous  for  its  diffraction  patterns,  and  so  on.  (2)  Cases  with  published  solutions  are 
preferred,  so  that  many  readers  will  already  be  familiar  with  the  problem.  The  intent 
is  to  maintain  the  focus  on  the  exact  ray  theory  developed  here,  and  demonstrations 
of  the  theory  are  much  more  comprehensible  when  applied  to  familiar  problems. 
(3)  Each  case  has  an  exact  analytic  solution  or  a  good  approximate  solution  suitable 
for  computation.  Numerical  solutions  of  the  Helmholtz  equation  are  avoided.  This 
gives  assurance  that  any  unexpected  behavior  of  the  rays  is  in  fact  valid  and  cannot 
be  attributed  to  the  vagaries  of  numerical  computation.  To  be  sure,  a  Runge-Kutta- 
Fehlberg[21]  numerical  integration  code  is  used  to  integrate  the  ray  path  equation, 
but  many  years  of  experience  with  numerical  ray  tracing  have  proven  this  procedure 
to  be  very  reliable.[22]  -[24]  Acoustic  propagation  in  realistic  ocean  environments  is 
an  extremely  complicated  process;  sophisticated  numerical  computer  models  often 
must  be  used  to  obtain  reasonably  accurate  simulations.  To  bring  such  models  into 
play  in  all  their  complexity,  while  simultaneously  trying  to  establish  the  fundamental 
properties  of  exact  ray  theory,  would  be  an  intolerable  distraction  for  the  purposes  of 
this  dissertation.  Instead,  we  treat  only  cases  where  the  acoustic  field  has  reasonably 
simple  behavior.  This  does  not  mean,  however,  that  exact  ray  analysis  is  limited  to 
“pencil  and  paper”  problems.  Some  numerical  propagation  models  can  be  modified 
to  compute  the  acoustic  field  gradients  as  well  as  the  field  itself  and  to  save  them  for 
use  by  a  ray  trace  program.  Such  models  would  then  be  able  to  predict  not  only  the 
amplitude  of  the  acoustic  field  but  also  its  direction  of  propagation. 

There  are  two  stages  to  the  construction  of  the  exact  ray  diagrams  which  follow. 
First,  the  Helmholtz  equation  is  solved  and  the  gradients  of  the  acoustic  field  are 
computed,  all  by  analytical  methods.  Then  the  ray  path  equation  dr /da  =  Im  V'J/4' 
is  solved,  usually  by  numerical  integration.  Since  the  solution  to  the  Helmholtz 
equation  is  different  for  each  case  and  usually  rather  complicated,  the  reader  will 
find  that  somewhat  lengthy  and  perhaps  tedious  mathematical  preliminaries  precede 
the  generation  of  each  ray  diagram.  I  have  tried  to  minimize  the  distraction  of  these 
mathematical  excursions  by  simply  providing  a  reference  if  a  suitable  solution  has 
been  published.  There  would  be  little  point  in  reproducing  solutions  to  problems 
selected  in  part  because  their  solutions  appear  in  the  literature.  Otherwise,  I  have 
tried  to  provide  sufficient  details  of  the  development  to  allow  reconstruction  of  the 
result  without  belaboring  the  derivation.  For  the  most  part,  this  material  can  be 
skimmed  with  little  loss  of  continuity  or  comprehension. 


As  a  simple  demonstration  of  exact  ray  theory  analysis  we  confirm  that  the  exact 
ray  trajectories  are  straight  lines  for  the  case  of  a  point  source  in  an  unbounded 
homogeneous  medium  (the  classical  ray  paths  are  also  straight  lines,  of  course). 

The  solution  to  the  Helmholtz  equation  for  a  point  source  is 


^  =  Qe^'^yr, 


Mm 


where  Q  is  the  source  strength,  with  dimensions  of  volume  per  unit  time,  and  r  is  the 
radial  distance  from  the  source,  which  is  assumed  to  lie  at  the  origin.  Upon  taking 
the  gradient  of  9  and  using  Vr  =  r/r,  we  find  that 


The  ray  paths  are  solutions  to  Eq.  (4-1 ); 


But  K  =  k  hy  Eq.  (4-3),  so  the  ray  path  equation  is 

dr  _  r 
ds  r 

In  words,  the  ray  path  tangent  vector  is  the  unit  radial  vector  and  the  rays  are 
therefore  radial  lines. 


C.  THE  LINEAR  SOUND  SPEED  CASE 

Classical  ray  theory  provides  an  exact  solution  to  the  Helmholtz  equation  when 
the  sound  speed  is  a  constant.  We  now  examine  the  linear  sound  speed  profile  case 
c  =  az,  where  classical  ray  theory,  with  a  slight  modification,  still  yields  an  exact 
solution.  Most  textbooks  on  underwater  acoustics  show  that  the  rays  for  a  linear 
profile  are  circular  arcs  centered  on  the  surface, [2 ,7, 10)  as  shown  in  Fig.  4.1. 
Classical  ray  analysis  of  this  problem  produces  the  approximate  solution 

_  2(zzo)»/» 

R1R2  ’ 

=  2-tanh-‘ (4-4) 

for  the  acoustic  field  '9  =  Ae'*,  where  Zq  is  the  source  depth  and 
Ri  =r^  -i-{z  -  Zof,  RI  =  +  (z  +  Zq}^. 

Pekeris[25]  showed  that  the  exact  solution  results  from  modifying  the  phase  (4-4)  to 

♦_.  =  2^U„h-(|).  (4-5) 

where 

a 


Exact  ray  analysis  of  the  two  preceding  cases  produced  results  which  were  little 
different  from,  or  not  at  all  different  from,  the  classical  ray  theory  results.  If  this 
always  occurred  then  exact  ray  theory  would  be,  at  most,  a  perturbative  correction 
to  classical  ray  theory.  We  now  treat  a  more  typical  problem,  henceforth  called  the 
Kormilitsin  problem  after  the  original  investigator,  where  the  classical  and  exact  ray 
diagrams  are  very  different.  The  case  fc*  =  klz/zo,  depicted  in  Fig.  4.2  along  with 
the  corresponding  sound  speed  profile,  is  notable  in  that  the  classical  ray  diagram 


exhibits  a  single  smooth  convex  caustic  and  so  serves  as  a  canonical  problem  in  the 
analysis  of  the  field  near  a  caustic. 


(a)  (b) 

Figure  4.2  The  Kormilitsin  profile.  The  monotonic  decrease  in  sound  speed 

WITH  DEPTH  CAUSES  DOWNWARD  REFRACTION  OF  BOTH  CLASSICAL  AND  EXACT  RAY  PATHS, 
(a)  (b)  c=  co(zoA)*/^. 


Kormilitsin  [26]  obtained  an  exact  elementary  solution  to  the  Helmholtz  equation 
for  a  line  source  at  zq;  the  corresponding  solution  for  a  point  source  is[27] 

* = (^)  T'’"’  h  + 4'" + -  %)] 


where  ko  =  k{zo).  Yet,  for  our  purposes,  a  good  approximation  which  lends  itself 
readily  to  computation  is  preferable  to  a  computationally  resistant  exact  solution. 
Ludwig[28]  and  Kravtsov[29]  independently  developed  a  uniform  zisymptotic  approx¬ 
imation  specifically  to  compute  the  field  near  a  simple  caustic;  it  is  therefore  the 
method  of  choice  for  this  problem.  In  deference  to  the  fact  that  the  uniform  asymp¬ 
totic  expansion  is  an  approximation,  the  rays  generated  from  it  will  be  referred  to  as 
frequency  dependent  rays  rather  than  exact  rays. 

The  development  of  several  wave  theoretical  corrections  to  ray  theory  was  mo¬ 
tivated  by  a  desire  to  retain  the  conceptual  and  computational  advantages  of  ray 
theory  while  obtaining  improved  or  corrected  predictions  where  the  ray  theory  ap¬ 
proximation  fails.  Ludwig’s  uniform  asymptotic  approximation  is  an  example  of 
such  a  correction  formula.  Beam  displacements,  discussed  in  Sec.  IV. H,  are  another. 


These  correction  formulas  typically  make  explicit  use  of  claissical  eigenray  path  infor¬ 
mation  already  computed  by  ray  models.  Ludwig’s  formula  makes  use  of  clcissical  ray 
amplitudes  and  phases  and  is  easier  to  understand  and  interpret  if  the  ray  analysis 
is  already  completed.  The  introduction  of  Ludwig’s  formula  is  therefore  postponed 
until  after  the  classical  ray  analysis  of  the  Kormilitsin  profile  problem  which  follows. 

1.  Classical  ray  analysis  of  the  Kormilitsin  profile 


The  classical  ray  analysis  of  the  Kormilitsin  profile  begins  with  the  deter¬ 
mination  of  the  ray  trajectories.  Once  the  ray  paths  are  known  the  ray  theoretical 
acoustic  field  amplitude  and  phase  along  a  ray  are  readily  found.  It  then  remains 
to  determine  the  rays  which  connect  source  and  observation  point;  these  are  the 
eigenrays. 

a.  Ra.y  trajectories 

We  obtain  the  classical  ray  paths  in  a  range  invariant  medium  by  convert¬ 
ing  Snell’s  law, 

cos*  0  ~  cos*  $0, 

into  a  differential  equation  whose  solution  reduces  to  quadrature.  Here,  6  is  the  ray 
path  angle  measured  with  respect  to  the  horizontal.  Snell’s  law  is  rendered  as 


cot^  = 


ko  cos  9o 

(F(2)- jtgcos^^o)*/* 


(4-6) 


with  the  aid  of  trigonometric  identities.  Since  9  is  the  ray  path  angle  {9o  is  the  path 
angle  at  the  source,  or  “launch  angle”),  the  derivative  of  the  ray  path  range  with 
respect  to  depth  is  given  by  drjdz  =  cot  9.  Equation  (4-6)  may  therefore  be  regarded 
as  a  first  order  differential  equation.  Moreover,  the  right  hand  side  of  Eq.  (4-6)  is  a 
function  of  z  only;  hence  Eq.  (4-6)  may  be  integrated  at  once  to  yield 


{k\z>)~klcos^eoyl'^ 

Setting  =  az  in  Eq.  (4-7)  and  integrating,  we  find 


r  =  2zy*(2*^*  sin0  —  2y*sin<?o)cos0 


Inverted,  this  is 


z  =  zo  +  r  tan  +  r*/ (4zo  cos*  Og). 


Thus  the  ray  paths  shown  in  Fig.  4.3  are  parabolic  arcs. 


For  a  point  source  on  the  z  axis  in  a  rcinge  invariant  medium,  the  classical 
ray  amplitude  is  obtained  from 

=  l/r|Cl,  (4-9) 

where  C  is  the  derivative  of  the  ray  path  depth  with  respect  to  launch  angle  at 
constant  range.  For  the  Kormilitsin  profile  this  is 


c  =  — 

^  aOo  r 
=  r  sec’  ( 


(4-10) 


From  Eq.  (4-10)  and  the  ray  path  equation  (4-8),  the  surface  C  =  0  is 


found  to  be 


r\  =  4zbZc. 


The  locus  of  points  (r,.,  defines  a  parabolic  surface  which  forms  an  envelope  to  the 
rays  and  separates  the  region  containing  rays  (the  insonified  region)  from  the  region 
devoid  of  rays  (the  shadow  zone),  as  shown  in  Fig.  4.3.  Such  surfaces  are  called 
caustics  because  of  the  unusually  intense  acoustic  fields  which  develop  along  them  at 
high  frequencies.  The  classical  ray  amplitude  is  singular  on  a  caustic. 


Figure  4.3  Classical  ray  diagram  for  the  Kormilitsin  profile.  The  parabolic  ray 
PATHS  (solid  curves)  ALL  EVENTUALLY  GRAZE  THE  CAUSTIC  (DASHED  CURVE),  WHICH  IS  ALSO 
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c.  Phase 


The  phase  in  the  ray  theory  approximation  is  given  by  the  path  integral 

$  =  r  kds\ 

Jo 


or 


$  =  /  k[l  +  {dz/dr'^Y^^  dr' 

Jo 

in  cylindrical  coordinates.  Using  k^  =  az  and  Eq.  (4-6),  the  phase  integral  becomes 

1/2 

1 

20> 


$=  ^  sec^o J  zir')dr'. 


Integration  of  the  ray  path  equation  (4-8)  then  gives 


$  =  k^r  sec 


l  +  ^tan9.  +  i(^)'sec’«„ 


(4-11) 


d.  Eigenrays 


Since  sec^  Aq  =  1  +  tan^  Bq,  Eq.  (4-8)  is  seen  to  be  quadratic  in  tan  Oq  as 
well  as  r.  The  two  roots  of  Eq.  (4-8)  for  a  point  (r,  z)  yield  the  launch  angles  of  the 
two  eigenrays,  or  rays  which  connect  the  source  at  (^o,^o)  to  the  point  (r,  2): 


(4-12) 


The  eigenray  launched  at  0q  grazes  the  caustic  before  reaching  the  observation  point, 
while  the  eigenray  launched  at  Oq  passes  through  the  observation  point  before  grazing 
the  caustic.  Henceforth,  all  ray  path  parameters  pertaining  to  the  eigenray  which 
gr«izes  the  caustic  before  reaching  the  observation  point  will  bear  the  superscript  1, 
while  path  parameters  of  the  second  eigenray  will  carry  the  subscript  2. 

Thus  the  acoustic  field  in  the  insonified  region  is  given  in  the  ray  theory 
approximation  by  the  coherent  sum  of  the  pressure  contributions  of  the  two  eigenrays: 

+  A2e'*\ 

where  the  amplitudes  A\  are  given  by  Eqs.  (4-9)  and  (4-10)  and  the  phases  $2  ^re 
given  by  Eq.  (4-11).  The  launch  angles  of  the  eigenrays  are  obtained  from  Eq.  (4-12). 

2.  Ludwig’s  uniform  asymptotic  approximation 

Having  obtained  the  classical  ray  theoretical  solution  to  the  Kormilitsin  prob¬ 
lem,  we  are  now  prepared  to  introduce  Ludwig’s  formula.  At  lowest  order  in  l/w, 
Ludwig’s  approximate  solution  to  the  Helmholtz  equation  is 

^  =  27r‘/^e-’^e*^[</+(r)  Ai(-e(r))  ig-{r)  Ai'(-^(r))], 


(4-13) 
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where 

i  =  l|(^‘  -  e  =  K**  +  *2),  fir*  =  ±  A,).  (4-14) 

Notice  the  use  of  the  classical  eigenray  amplitudes  and  phases  in  Exjs.  (4-14). 

The  Airy  function  Ai  is  a  solution  to  Airy’s  differential  equation [30] 

Ai"(i)  =  X  Ai(x). 

It  is  oscillatory  for  i  <  0  and  decays  exponentially  for  i  >  0. 

The  gradient  of  the  acoustic  field  is  needed  in  order  to  trace  frequency  de¬ 
pendent  rays.  After  calculating  the  gradients  of  Eiqs.  (4-13)  and  (4-14),  we  get 

Ai(-ori4-Ai^(-Ora 
^  g+  A\{-0  +  i9~  AiV-O’ 

where 

r,  =  ig^Ve-\-Vg^  +  i(g-V^, 
r2  =  -g~ve +  iVg~  -  g'^Vi, 

and 

-  ^^2)> 

V0  =  |(V$‘  +  V$2), 

Before  computing  the  gradients  of  the  cl2issical  ray  quantities  appearing  in 
Eqs.  (4-16),  it  proves  convenient  to  define  the  following  variables: 

I  =  z/zq,  m  =  r/2zQ,  =  (1  -b  m  tan  ^o)^- 

Then  the  ray  path  equation  may  be  written  compactly  as 

I  =  rn^  +  n*, 

the  eigenray  amplitudes  are 

*  2zo[n(l  +  /±2n)]>/2’ 

and  the  phases  are 

$‘  =  |ifco^o(l  +  fT^)(l  +  /±2n)*'^ 

The  gradients  of  the  eigenray  amplitudes  and  phases  appearing  in  Eq.  (4-16)  take 
the  forms 

VA\  =  -(A‘)^[e,m(l  -|-/±4n)  -e,(l  -|-/±4n  +  2n2)]  (4-17) 

n 

and 


(4-15) 


(4-16) 
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3.  Frequency  dependent  ray  analysis  of  a  caustic 

With  Eqs.  (4-15)  -  (4-18)  in  hand  for  a  frequency  dependent  ray  dia¬ 

gram  is  constructed  by  numerically  solving  the  ray  path  equation  drldcr  =  Im(  V'l'/^') 
for  various  launch  angles,  using  mathematical  software  library  routines  to  compute  Ai 
and  Ai'.  These  routines  must  be  capable  of  computing  Airy  functions  of  complex  ar¬ 
gument  if  the  shadow  zone  is  to  be  treated. 

Figure  4.4  shows  the  frequency  dependent  ray  diagrams  for  source  frequen¬ 
cies  of  5  Hz  and  50  Hz,  with  the  source  located  1000  m  below  the  surface  and  a  sound 
speed  at  the  source  of  1500  m/s.  These  frequency  dependent  rays  obviously  differ 
starkly  from  their  classical  counterparts  shown  in  Fig.  4.3.  The  frequency  dependent 
ray  diagrams  can  be  divided  roughly  into  four  zones  according  to  the  behavior  of  the 
rays  found  within  them.  Near  the  source  the  rays  form  the  radial  lines  characteristic 
of  spherical  spreading,  as  discussed  in  Sec.  IV.B.  In  the  insonified  region  far  from  the 
source  but  well  inside  the  caustic,  the  frequency  dependent  rays  begin  to  undulate. 
As  the  rays  approach  the  caustic,  the  undulations  damp  out  and  the  rays  bunch  to¬ 
gether,  running  parallel  to  each  other  and  almost  parallel  to  the  caustic.  The  greatest 
constriction  of  the  ray  bundles  is  found  somewhat  inside  the  caustic.  Eventually,  the 
rays  cross  the  caustic  and  intrude  into  the  shadow  zone,  where  they  begin  to  spread 
out.  And,  in  sharp  contrast  to  the  classical  rays,  these  rays  never  cross.  In  fact  they 
fill  all  of  space,  so  that  for  every  observation  point,  even  in  the  shadow  zone,  there 
is  exactly  one  eigenray. 

The  details  of  these  strange  behaviors  are  sensitive  to  the  source  frequency. 
The  amplitudes  of  the  undulations  diminish  with  increasing  frequency,  but  the  un¬ 
dulations  become  more  rapid,  giving  the  50  Hz  rays  a  “jittery”  appearance.  Even 
though  they  were  launched  at  the  same  angles  as  the  5  Hz  rays  (from  —85°  to  -|-85° 
in  5°  increments),  the  50  Hz  rays  are  much  more  concentrated  near  the  caustic  and 
reach  their  greatest  concentration  nearer  the  caustic.  Nor  do  they  venture  nearly  as 
far  out  into  the  shadow  zone  as  the  5  Hz  rays. 

The  frequency  dependent  ray  diagrams  can  be  used  to  explain  the  behav¬ 
ior  of  propagation  loss  curves.  Horizontal  reference  lines  are  drawn  across  both 
ray  diagrams  in  Fig.  4.4  at  the  1200  m  depth  level.  Figure  4.5  shows  the  propa¬ 
gation  loss  curves  at  this  depth  for  both  frequencies,  computed  by  Ludwig’s  uni¬ 
form  asymptotic  approximation  (solid  curves)  and  by  claissical  ray  theory  (dashed 
curves).  (The  propagation  loss  curves  plotted  here  were  computed  using  the  for¬ 
mula  PL  =  — 10  logjo[^^*].)  The  ray  theory  propagation  loss  curve  is  extended  into 
the  shadow  zone  by  a  slight  modification  of  classical  ray  theory  to  be  described  later. 
The  ray  theory  results  agree  very  closely  with  the  uniform  asymptotic  approximation 
except  in  the  vicinity  of  the  caustic. 

The  separation  between  the  rays  as  they  cross  the  reference  .ine  is  a  measure  of 
the  spreading  of  the  rays  and,  hence,  the  intensity  of  the  acoustic  field.  Comparing 
the  5  Hz  ray  diagram  with  the  corresponding  propagation  loss  curve,  we  see  that 
wherever  the  undulations  of  a  pair  of  rays  causes  them  to  draw  apart  as  they  cross 
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Figure  4.4  Frequency  dependent  ray  diagram  for  the  Kormilitsin  profile,  (a)  At 
5  Hz  THE  RAYS  undulate  GENTLY  AND  READILY  INTRUDE  INTO  THE  SHADOW  ZONE,  (b)  AT 
50  Hz  THE  UNDULATIONS  ARE  MUCH  MORE  RAPID  AND  OF  SMALLER  AMPLITUDE.  PENETRATION 


OF  THE  SHADOW  ZONE  IS  GREATLY  REDUCED. 
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Figure  4.5  Propagation  loss  curves  for  the  Kormilitsin  profile.  The  broken  verti¬ 
cal  LINE  INDICATES  THE  RANGE  OF  THE  CAUSTIC  AT  Z  =  1200  m.  ThE  SOURCE  IS  AT  1000  m, 
OBSERVATION  POINTS  ARE  AT  1200m.  (a)  5Hz.  (b)  50 Hz. 
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the  reference  line  there  is  a  dip  in  the  field  intensity  over  the  corresponding  range 
interval;  similarly,  the  intensity  is  greatest  where  the  ray  pair  draws  together.  The 
undulations  vanish  as  the  reference  line  nears  the  caustic;  the  rays  become  most 
concentrated  just  inside  the  caustic  and  then  begin  to  spread  out  as  the  reference 
line  passes  through  the  caustic  and  into  the  shadow  zone.  Likewise,  the  intensity 
reaches  its  maximum  at  the  point  of  maximum  concentration  of  the  rays  and  then 
decays  monotonically  as  the  range  increases  to  and  beyond  the  range  of  the  caustic. 

The  50  Hz  ray  diagrams  and  propagation  loss  curves  possess  the  same  general 
features  as  those  for  5  Hz  but  differ  in  the  details.  The  rapid  undulations  of  the  rays 
correspond  to  the  rapid  oscillations  in  the  field  intensity.  The  heavier  concentration 
of  rays  nearer  the  caustic  and  the  lesser  degree  of  intrusion  into  the  shadow  zone  co¬ 
incides  with  the  peak  intensity  occurring  nearer  the  caustic  and  the  rapid  exponential 
decay  of  the  field  in  the  shadow  zone. 

The  propagation  loss  curves  are  interesting  in  their  own  right  and  deserve  a 
few  remarks.  Comparison  of  ray  theory  with  uniform  asymptotic  predictions  confirms 
the  maxim  that  ray  theory  remains  valid  in  the  insonified  region  up  to  the  last  zone  of 
constructive  interference  just  inside  the  caustic.  The  acoustic  field  intensity  reaches 
a  local  maximum  in  this  zone  and  then  begins  to  decay.  The  intensity  on  the  caustic 
itself  is  actually  several  decibels  below  the  local  maximum. 

But  notice  that  the  “classical”  ray  theory  propagation  loss  curve  continues 
into  the  shadow  zone  even  though  the  classical  ray  diagram  of  Fig.  4.3  clearly  shows 
the  rays  confined  to  the  insonified  region.  This  portion  of  the  propagation  loss  curve 
was  generated  using  the  complex  ray  theory  of  White  and  Pedersen. [31]  We  had 
assumed  until  now  that  only  real-valued  eigenray  launch  angles  and  ray  path  coor¬ 
dinates  are  admissible.  But  the  ray  path  equation  (4-8)  is  quadratic  in  tan  do  S’nd 
yields  complex- valued  eigenray  launch  angles  for  observation  points  in  the  shadow 
zone.  The  path  parameters  of  rays  with  complex  launch  angles  take  on  complex 
values.  Ludwig’s  asymptotic  formula  (4-13)  remains  valid  in  the  shadow  zone  and 
makes  explicit  use  of  the  complex  amplitudes  and  phases  of  the  complex  eigenrays. 
Following  White  and  Pedersen,  we  compute  the  field  in  the  shadow  zone  by  retaining 
the  complex  eigenray  labeled  2,  but  discard  the  eigenray  labeled  1  because  it  gives 
rise  to  a  nonphysical  exponentially  growing  field.  In  other  words  =  A^e'*'  +A2€'*^ 
in  the  insonified  region,  but  ^  =  /I2C***  in  the  sheulow  zone. 


Unlike  the  first  two  cases  examined  in  this  chapter,  the  frequency  dependent 
rays  for  the  Kormilitsin  profile  are  irreconcilably  different  from  their  classical  coun¬ 
terparts;  yet  both  kinds  of  rays  are  purported  to  represent  the  paths  of  energy  propa¬ 
gation.  It  seems  appropriate  to  try  to  resolve  the  paradox  before  proceeding  to  more 
demonstrations  of  exact  ray  theory. 

The  solution  to  the  dilemma  lies  in  the  nature  of  classical  ray  theory  and  the 
ray  theory  approximation.  Classical  ray  theory  has  been  used  for  many  d<;cades  to 
study  wave  propagation  and  researchers  have  thereby  gained  a  great  deal  of  expe- 
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rience  and  intuition  which  guides  their  expectations  of  the  behavior  of  rays.  It  is 
therefore  disconcerting  to  many  people  to  be  confronted  with  exact  ray  diagrams 
which  systematically  violate  their  expectations  of  what  ray  diagrams  should  look 
like.  It  is  tempting  to  dismiss  these  curves  as  being  something  other  than  rays  (and, 
to  be  sure,  they  are  indeed  not  classical  rays).  But,  although  it  is  certainly  correct 
to  claim,  for  example,  that  these  paths  are  the  orthogonal  trajectories  of  the  solution 
to  the  Helmholtz  equation,  this  assertion  does  not  serve  to  distinguish  these  paths  as 
something  other  than  rays.  These  are  precisely  the  paths  which  result  from  avoiding 
the  approximation  invoked  in  the  Sommerfeld  derivation  of  the  classical  ray  path 
equation;  they  form,  now  with  perfect  fidelity,  the  energy  flux  tubes  which  are  the 
essential  feature  of  a  ray  theory.  In  short,  these  paths  are  rays. 

Another  tempting  resolution  of  the  par2Miox,  also  incorrect,  is  to  suppose  that 
if  exact  rays  were  constructed  at  high  enough  frequencies  then  they  would  come 
to  resemble  classical  rays,  because  claissical  ray  theory  is  a  high  frequency  approxi¬ 
mation.  But  this  expectation  will  usually  be  frustrated,  if  for  no  other  reason  than 
because  exact  ray  diagrams  never  exhibit  multipathing  no  matter  how  high  the  source 
frequency,  whereas  multipathing  is  a  characteristic  feature  of  most  classical  ray  di¬ 
agrams.  Classical  ray  theory  is  indeed  a  high  frequency  asymptotic  approximation, 
but  the  asymptotic  expansion  is  nonuniform  in  most  cases,  and  the  nonuniformity 
manifests  in  the  vicinity  of  caustics.  Caustics  are  an  inevitable  consequence  of  mul¬ 
tipathing;  if  the  ray  paths  cross  at  all  then  there  will  usually  be  a  limiting  surface 
where  the  rays  cross  tangentially.  Such  a  surface  is,  of  course,  a  caustic,  and  the  dif¬ 
ferential  ray  tube  cross  sectional  areas  vanish  on  the  caustic.  The  transport  equation 
then  predicts  nonphysical  singular  acoustic  fields  on  the  caustic. 

When  the  classical  ray  path  equation  is  derived  by  the  Sommerfeld  development 
rather  than  as  a  high  frequency  approximation,  one  makes  the  ray  theory  assumption 
that  V*A//4  <C  But  this  assumption  usually  turns  out  to  be  much  more  than  a 
perturbative  approximation;  the  invocation  of  the  ray  theory  approximation  changes 
the  resulting  rays  from  exact  rays  to  classical  rays  which  may  behave  quite  differently. 

Yet  classical  ray  theory  has  been  used  to  solve  the  Helmholtz  equation  for  many 
years  and  with  considerable  success.  If  further  proof  were  needed,  the  propagation 
loss  curves  computed  for  the  Kormilitsin  profile  using  classical  ray  theory  agree  ex¬ 
tremely  well  with  those  computed  by  Ludwig’s  uniform  asymptotic  approximation, 
except  near  the  caustic.  Classical  ray  theory  was  clearly  not  “wrong”  in  that  case, 
so  why  were  the  classical  rays  so  uncompromisingly  different  from  the  frequency  de¬ 
pendent  rays?  In  particular,  why  is  the  phenomenon  of  multipathing,  so  familiar  in 
classical  ray  diagrams,  absent  from  exact  ray  diagrams? 

In  order  to  solve  the  multipathing  conundrum,  we  revisit  the  Kormilitsin  prob¬ 
lem,  but  now  we  consider  a  quiescent  medium  disturbed  by  a  point  source  which  is 
suddenly  activated.  The  reader  will  recall  from  the  discussion  in  Sec.  III.D.l  that 
the  leading  edge  of  the  resulting  disturbance  propagates  according  to  the  predictions 
of  classical  ray  theory.  Specifically,  a  listener  in  the  insonified  region  will  note  two 
distinct  arrivals  at  times  and  intensities  predicted  by  ray  theory,  the  first  arrival  be- 


ing  from  that  portion  of  the  surface  formed  by  the  propagating  discontinuity  which 
encountered  the  caustic  before  sweeping  past  the  listener.  The  predictions  of  classical 
ray  theory  are  also  in  accord  with  common  experience  of  multipathing  in  the  form  of 
echoes.  Ray  theory  is  unsurpassed  at  tracking  propagating  discontinuities. 

Multipathing,  whether  due  to  reflection  or  refraction,  is  a  quintessentially  broad¬ 
band  phenomenon.  But  a  listener  exposed  to  a  continuous  wave  (cw)  source  hears 
only  a  pure  tone;  there  is  no  sensation  of  distinct  arrivals — no  multipathing.  Whereas 
classical  ray  theory  describes  the  broadband  phenomenology  to  which  it  is  better 
suited,  exact  ray  theory  describes  cw  phenomenology,  and  the  two  are  simply  differ¬ 
ent. 

One  can  reconstruct  the  acoustic  disturbance  created  by  a  broadband  source  by 
Fourier  synthesis  of  the  acoustic  flelds  generated  by  cw  sources,  each  radiating  at 
a  different  frequency.  The  acoustic  flelds  of  each  cw  source  are  best  represented  by 
exact  rays,  while  the  broadband  disturbance  is  best  represented  by  classical  rays. 
When  one  applies  classical  ray  theory  to  the  Helmholtz  equation,  one  presses  into 
service  the  solution  to  a  broadband  problem  in  order  to  solve  a  cw  problem.  This  is 
permissible — at  the  cost  of  failures  of  the  theory  near  caustics,  as  is  well  understood. 
It  is  less  well  understood  that  classical  ray  theory  comes  freighted  with  the  conceptual 
baggage  of  broadband  propagation.  Ironically,  much  of  the  conceptual  framework 
and  many  intuitive  notions  of  the  behavior  of  cw  fields  are  based  on  experience  with 
classical  ray  theory. 

F.  PROPAGATION  IN  A  WAVEfltJlDR 

Next  we  treat  the  case  where  k^{z)  is  taken  to  be  the  quadratic  function 

k‘^  =  kl[l-il-z/z,Y].  (4-19) 

Thus  k  has  a  maximum  value  of  ka  z  =  and  vanishes  at  «  =  0  and  z  =  2z„,  as 
shown  in  Fig.  4.6(a).  Figure  4.6(b)  shows  that  the  corresponding  sound  speed  profile 
c{z)  =  u:lk{z)  has  a  minimum  value  of  Ca  at  the  axis  depth  Za  and  increases  without 
bound  near  the  asymptotes  at  z  =  0  and  z  =  2za.  Sound  tends  to  refract  away  from 
regions  of  high  sound  speed  and  toward  regions  of  low  speed.  This  medium  therefore 
acts  as  a  waveguide  even  though  it  lacks  true  boundaries,  turning  both  classical  and 
exact  rays  away  from  the  asymptotes  and  steering  the  acoustic  energy  towards  the 
sound  channel  axis.  The  waveguide  thus  formed  turns  out  to  be  the  acoustic  analog  of 
the  quantum  harmonic  oscillator  and  henceforth  will  be  called  the  harmonic  oscillator 
waveguide. 


The  classical  ray  paths  may  be  determined,  as  in  the  case  of  the  Kormilitsin 
problem,  by  integrating  Snell’s  formula: 


Figure  4.6  The  harmonic  oscillator  waveguide.  (a)  it- 
(b)  c  =  Ca[l  -  (1  - 


k-l[l  -  (1  - 


Quadrature  is  simplified  when  the  source  is  located  on  the  axis.  Not  surprisingly,  in 
view  of  the  nature  of  the  waveguide,  the  ray  paths  turn  out  to  be  sinusoidal: 


r  =  To  +  x^a  sin  sin 


The  rate  of  change  of  ray  path  depth  with  respect  to  launch  angle  is 


(>  = 


dz 


=  Xacos^osinf - p  )  +  rtan^^ocos  ( - — :r  )  • 

\  Xa  cos  BqJ  \Za  COS  9q  ) 

Caustics  are  located  where  (  =  0.  For  a  given  launch  angle  Oq.  the  ranges  r.(0o)  at 
which  the  ray  grazes  the  caustic  turn  out  to  satisfy  the  nonlinear  equation 


tan(rc/;a  cos^o) 


=  —  tan^  00- 


rjz,  cos9o 

while  the  corresponding  depths  Xc(rc.^o)  satisfy  the  ray  path  equation 

Xc  =  Xo  +  Xa  sin  $0  sin  (  — )  . 

\  Xa  cos  9o  / 

2.  Normal  mode  representation  of  the  field  in  a  waveguide 


Waveguide  problems  can  often  be  solved  conveniently  using  normal  modes  (a 
brief  review  of  the  method  follows),  although  the  normal  modes  and  eigenvalues  must 
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usually  be  determined  numerically.  The  profile  (4-19)  selected  for  this  example  is  one 
of  very  few  problems  in  which  the  normal  modes  and  eigenvalues  can  be  expressed 
analytically  and  computed  easily. 

The  Green  function  ro)  satisfies  the  inhomogeneous  Helmholtz  equation 
V^'Pg  +  =  -^(r  -  ro). 


In  cylindrical  coordinates  this  is 


I  d  f  d\  a"  1  1  T  1.2  r  1  r/ 


ro)6{z-zo)£{<i>-<f>o).  (4-20) 


Because  the  source  is  located  on  the  z  axis  and  the  medium  is  horizontally  strati¬ 
fied,  the  acoustic  field  is  azimuthally  symmetrical  and  derivatives  with  respect  to  <p 
consequently  vanish.  The  inhomogeneous  Helmholtz  equation  (4-20)  becomes 


+  t*a-  --6{r)S{z  -  Zo)SW. 


(4-21) 


Integration  of  Eq.  (4-21)  with  respect  to  ^  from  0  to  27r  produces  the  desired  form 
of  the  inhomogeneous  Helmholtz  equation: 


d^'^G  Id^G  d^^G 

dr"^  r  5r  **"  dz"^ 


+  +  ^  +  k^^G  =  -~6{r)6{z  -  zo). 


(4-22) 


Since  the  homogeneous  Helmholtz  equation  corresponding  to  Eq.  (4-22)  is 
separable,  it  is  convenient  to  express  the  solution  to  the  inhomogeneous  Helmholtz 
equation  in  terms  of  the  solutions  to  the  separated  depth  and  radial  equations.  It 
can  be  shown  that  the  Green  function  for  Eq.  (4-22)  is  given  by  (see  Boyles, [2] 
pp.  163-177) 


^g(»',  Z;  2o)  =  J  XI  «n(zb)«n(2)  H^^\knr). 

**  n=0 


(4-23) 


The  radial  function  Ho\  the  function  of  the  first  kind,  is  a  solution  of  a 
Bessel  equation  which  may  be  written  in  the  form 

d 

r^—Rc{Kr)  4-  r— i?o('«»’)  -|-  «*r*ilo(«:r)  =  0. 
ar'  ar 

//o*^(«:r)  is  that  solution  to  the  Bessel  equation  which  is  singular  at  r  =  0  and  which 
is  consistent  with  outgoing  waves  and  an  acoustic  field  with  time  dependence  e"*"'; 
that  is,  it  satisfies  the  Sommerfeld  radiation  condition 


lim 


(ml 

^  dr 


■(«r)  -  ik„  H^^\Kr) )  =  0 


ft 
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The  normal  modes  u„(ir)  are  the  eigenfunctions  of  the  depth-separated  Helm¬ 
holtz  equation 

^  +  (t*  -  kl)u„  =  0,  (4-24) 

and  the  awe  the  corresponding  eigenvalues.  The  normal  modes  satisfy  the  or¬ 
thonormality  condition 

/oo 

u„,{z)un{z)  dz  =  (5„„, 

•OO 

and  boundary  or  radiation  conditions  imposed  on 't  and  d^fdz  at  horizontal  bound¬ 
ary  surfaces. 

In  the  case  of  the  harmonic  oscillator  waveguide,  the  normal  modes  are  ex¬ 
pressible  in  terms  of  Hermite  polynomials[30]: 

Un{z)  = 

^  -  Za), 

and  the  eigenvalues  are  equally  spaced  in  increments  of  2kalza‘. 

kl^kl-{2n-^l)kjz,.  (4-25) 

The  Hermite  polynomials  Hn  (not  to  be  confused  with  Hankel  functions,  such  as 
satisfy  the  differential  equation 

//:'(x)-2xi/;(x)-|-2n7/„(x)=0, 

and  the  normalization  condition 

r  e-^[Hr,ix)]Ux  =  2’*T*/*n! . 

J—OO 

They  are  easily  computed,  along  with  their  derivatives,  using  the  recursion  relations 

Hn{x)  =  2xHn-i{x)  -  2(n  -  l)Hn-2{x), 

Hoix)  =  1, 

Hi{x)  =  2x, 

and  the  differentiation  formula 


/f'(x)  =  2nHn-\{x).  (4-26) 

One  of  the  properties  of  a  waveguide  is  that  not  all  modes  propagate  freely 
within  it;  beyond  some  mode  cutoff  number  Me  the  modes  are  damped  exponentially 
as  they  travel  down  the  waveguide.  We  exploit  this  property  to  terminate  the  normal 
mode  summation  (4-23).  For  the  harmonic  oscillator  waveguide,  mode  cutoff  occurs 
at 


"c  =  -  1)J- 


(4-27) 
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The  eigenvalues  are  negative  for  modes  with  n  >  ric,  as  one  may  confirm  by  inspection 
of  the  dispersion  relation  (4-25).  The  argument  of  the  Hankel  function,  is 
consequently  imaginary  for  modes  above  cutoff.  For  large  argument  w,  the  Hankel 
function  has  the  asymptotic  form 

V  vw 

hence,  for  a  purely  imaginary  argument  w  =  ix. 


Modes  above  cutoff  are  thus  exponentially  damped  and  do  not  propagate  very  far 
down  the  waveguide.  One  may  therefore  truncate  the  mode  summation  at  cutoff, 
provided  the  calculations  are  restricted  to  ranges  beyond  the  effective  propagation 
distance  of  the  first  damped  mode. 

We  are  not  quite  ready  to  trace  rays;  we  still  need  V'J.  The  gradient  of  the 
mode  sum  (4-23)  for  'i'  is 

0 

•  oo 

=  Ho\knr)  -  er*:„u„(z) //j’^(*:„r)]. 

(4-28) 

Using  Eq.  (4-26)  for  the  derivative  of  a  Hermite  polynomial,  we  obtain 

«n(^)  =  -  xHn{x)].  (4-29) 

With  Eq.  (4-23)  for  'I'g,  and  Eqs.  (4-28)  and  (4-29)  for  V^g,  the  ray  trace  proceeds 
by  solving  dr/</<7  =  Im(V^/$)  numerically. 

3.  Exact  ray  analysis  of  a  refractive  waveguide 

In  this  example  we  take  the  axis  depth  Za  to  be  500  m  and  the  sound  speed 
on  the  axis,  to  be  1500  m/s.  For  a  source  frequency  of  3  Hz,  mode  cutoff  occurs 
at  Tic  =  L’"’  ~  V2J  =2  by  Eq.  (4-27);  the  n  =  2  mode  propagates  freely  while 
the  n  =  3  mode  attenuates  exponentially  with  range.  The  low  frequency  and  corre¬ 
spondingly  low  mode  cutoff  number  ensure  that  the  acoustic  field  will  have  simple 
behavior,  and  that  any  instabilities  in  the  Hermite  polynomial  recursion  relation  will 
not  have  an  opportunity  to  cause  an  accumulation  of  errors  that  might  compromise 
the  accuracy  of  the  field  calculation. 

Figure  4.7  compares  classical  rays  with  exact  rays  for  a  source  located  on 
the  sound  channel  axis  [zq  =  Za  =  500  m).  Since  modes  above  cutoff  may  not  be 
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Figure  4.7  Exact  rays  and  classical  rays  in  the  harmonic  oscillator  waveguide. 
(a)  The  classical  rays  (solid  curves)  form  an  infinite  series  of  cusped  caustics 
(dashed  curves),  (b)  The  exact  rays  never  cross  and  therefore  never  form  caustics. 
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dropped  from  the  mode  summation  until  they  have  traveled  far  enough  down  the 
waveguide  to  damp  out,  exact  ray  traces  must  be  confined  to  regions  far  enough 
from  the  source  that  the  mode  sunrunation  can  be  conveniently  truncated  with  an 
acceptable  truncation  error.  The  exact  ray  traces  depicted  in  Fig.  4.7(b)  commence 
1  km  downrange  from  the  source.  The  launch  angles  of  these  rays  are  not  known 
and  cannot  readily  be  identified  with  their  classical  counterparts  of  Fig.  4.7(a);  the 
ability  to  initiate  the  exact  ray  trace  away  from  the  source,  as  discussed  in  Sec.  IV. A, 
is  indispensible  here. 

The  classical  and  exact  rays  are  once  again  uncompromisingly  different.  Just 
as  in  the  Kormilitsin  problem  of  Sec.  IV.D,  the  exact  rays  do  not  cross  and  form 
caustics,  although  they  do  concentrate  in  regions  of  high  intensity.  The  classical  rays 
do  cross  systematically  and  thus  form  an  infinite  series  of  cusped  caustics,  as  shown 
in  Fig.  4.7(a).  The  regions  of  greatest  constriction  of  the  exact  ray  bundles  lie  on 
the  sound  channel  aocis,  but  do  not  always  coincide  with  the  cusps  of  the  classical  ray 
caustics,  where  the  acoustic  field  is  most  intense  at  high  frequencies. 

4.  Phase  velocity 

The  complex  velocity  potential  is  explicitly  time  dependent. 

The  phase  speed  Cph  is  defined  to  be  the  rate  ds/dt  at  which  a  point  must  advance 
along  a  ray  path  in  order  to  keep  the  phase  of  (p  constant  at  that  point,  where  s  is 
distance  along  the  ray  path.  By  requiring  that  the  time  derivative  of  ^[s(<)]  —  ujt 
vanish  and  by  using  d^/ds  =  K,  the  phase  speed  is  determined  to  be 

Cph  =u;/K. 

In  the  ray  theory  approximation  fc  w  /C, 

Cph  «  u?/fc  =  c. 

In  view  of  Eq.  (4-30),  the  exact  ray  path  equation 


may  be  rewritten  as 

±(±±]=v(-L 

ds  y  Cph  ds  J  y  Cph 

Conventional  ray  trace  programs  are  based  on  the  assumption  that  the  sound  speed 
and  phase  speed  are  synonymous;  they  solve  Eq.  (4-31)  by  assuming  that  Cpi,  =  c. 
But  even  a  conventional  ray  trace  program  would  nevertheless  trace  out  the  exact 
ray  paths  if  it  were  supplied  the  true  phase  speed  (4-30)  in  place  of  the  sound  speed. 

Figure  4.8  shows  the  level  contours  of  the  phase  speed  in  the  harmonic  oscil¬ 
lator  waveguide  for  different  source  depths  and  frequencies.  The  phase  speeds  were 
calculated  using  Cph  =  with  I\  obtained  from  Eq.  (4-3).  Phase  speed  contours 


(4-31) 


(4-30) 


48 


for  the  3  Hz  source  at  500m  are  shown  in  Fig.  4.8(a).  Notice  that  the  phase  speed 
is  markedly  range  dependent  even  though  the  sound  speed  is  range  invariant.  The 
phase  speed  contours  shift  when  the  frequency  is  raised  slightly  from  3.0  Hz  to  3.1  Hz 
(Fig.  4.8[b]),  and  when  the  source  depth  is  moved  down  slightly  from  500  m  to  525  m 
(Fig.  4.8[c]). 

In  general,  K  and  the  phase  speed  also  depend  on  the  source  configuration 
and  on  the  boundaries  or  radiation  conditions.  The  frequency  dependence  of  the 
phase  speed  is  concomitant  with  the  frequency  dependence  of  the  exact  ray  paths 
and  with  the  dispersion  observed  in  inhomogeneous  media,  even  in  the  absence  of 
boundaries. 

Wc  conclude  the  discussion  of  the  harmonic  oscillator  waveguide  with  a  fur¬ 
ther  remark  on  the  origin  of  its  name.  Hermite  polynomials  appear  in  eigenfunction 
solutions  to  the  Schrodinger  equation  describing  a  quantum  particle  trapjied  in  a 
conservative  potential  where  the  restoring  force  on  the  particle  is  proportional  to 
its  displacement  from  equilibrium  (almost  all  quantum  mechanics  textbooks  discuss 
this  famous  problem).  In  the  classical  limit,  the  particle  undergoes  time  harmonic 
oscillations  about  its  equilibrium  position.  The  author  chose  the  name  “harmonic 
oscillator”  waveguide  to  describe  the  waveguide  which  is  the  acoustic  analog  of  the 
quantum  harmonic  oscillator. 

G.  FRAUNHOFER  DIFFRACTION 

We  now  take  up  the  problem  of  a  circular  piston  of  radius  a  mounted  flush  in  an 
infinite  plane  baffle.  The  farfield  acoustic  radiation  exhibits  Fraunhofer  diffraction 
patterns. 

1.  Radiation  from  a  circular  piston  in  a  baffle 

In  an  unbounded  homogeneous  medium,  the  acoustic  field  at  r,  due  to  a  time 
harmonic  point  source  at  Tq,  is 

_«Ar|r-ro| 

^G(r;ro)  =  g^ - :•  (4-32) 

ir-ro| 

To  a  good  approximation,  each  elementary  surface  on  the  piston  face  radiates  accord¬ 
ing  to  Eq.  (4-32),  so  that  the  total  acoustic  field  is  given  by  an  integration  over  the 
surface  of  the  piston  (we  are  exploiting  the  Huygens- Fresnel  principle).  The  integra¬ 
tion  is  most  readily  performed  with  the  system  configuration  and  coordinate  system 
shown  in  Fig.  4.9.  The  baffle  and  the  surface  of  the  piston  lie  in  the  xy  plane  with 
the  center  of  the  piston  at  the  origin.  The  piston  radiates  upward  into  the  z  >  0 
region. 

The  integration  over  the  piston  surface  is  accomplished  by  resort  to  a  farfield 
approximation.  Tiie  diffraction  patterns  manifest  only  in  the  farfield,  where  r  a. 
The  details  of  the  calculation  are  given  in  many  references  on  optics  and  acoustics 
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Figure  4.9  The  piston  in  a  baffle.  A  circular  piston  mounted  in  a  plane  baffle 
CREATES  Fraunhofer  diffraction  patterns. 


(see  Kinsler  and  Frey[32]  for  an  elementary  discussion).  The  result  is 

2Ji{kas\nO) 


^'(r)  =  Qtto^ 


ka  sin  6 


(4-33) 


for  the  acoustic  field  far  from  the  piston. 

The  directivity  factor  2  Ji{ka  sin  B)/{ka  sin  0)  accounts  for  the  Fraunhofer 
diffraction  zones.  The  function  2  Ji(i)/i  is  oscillatory,  with  an  amplitude  that  falls 
off  as  for  large  x.  As  expected,  a  beam  of  sound  appears  directly  in  front  of 
the  piston,  where  0  =  0  and  2  Ji{ka  sin  9)/ {ka  sin  0)  attains  its  maximum  value  of 
unity.  But  if  the  piston  is  large  compared  to  a  wavelength  {ka  1),  then  conical 
zones  of  enhanced  intensity,  due  to  diffraction,  also  form  at  angles  corresponding  to 
the  extrema  of  2  Ji(fcasin^)/(/:asin0),  interspersed  with  shadow  zones  correspond¬ 
ing  to  the  zeros  of  Ji{kasin0).  On  the  other  hand,  if  the  piston  is  on  the  order  of  a 
wavelength  across  {ka  ~  1),  then  no  shadow  zones  or  secondary  maxima  appear.  If 
the  piston  is  small  {ka  ■Cl),  then  the  primary  beam  is  greatly  dispersed,  radiating 
with  nearly  uniform  intensity  in  all  directions  from  the  piston. 

When  secondary  maxima  are  present  they  are  far  less  intense  than  the  pri¬ 
mary.  At  any  fixed  distance  r  from  the  piston,  the  peak  intensity  in  the  first  diffrac¬ 
tion  zone  is  17.5dB  below  that  of  the  central  maximum;  the  peak  intensity  of  the 
second  zone  is  23.8  dB  down. 


2.  Diffraction  of  exact  rays 

Upon  taking  the  gradient  of  ^  (Eq.  [4-33]),  we  discover  that 


V'if 


which  is  the  same  result  obtained  in  the  case  of  the  point  source  in  an  unbounded 
hoinog^.»^ous  medium  (see  Sec.  IV. B).  The  reader  will  recall  that,  in  that  case,  the  ray 
paths  were  radial  lines  originating  at  the  origin.  Here,  our  conclusions  are  rt'stricted 
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to  the  farfield,  where  the  ray  paths  are  once  again  found  to  be  radial  lines,  with 
retrograde  projections  back  to  the  origin.  (However,  we  do  not  know  the  paths  the 
rays  followed  through  the  nearfield  to  get  to  the  farfield.) 

The  field  amplitude,  obtained  from  A?  =  4'^*,  is 

Ji(A:asin0) 
ka  sin  0)  ’ 

which  is  to  be  compared  with  A  =  l/r  for  the  point  source.  Clearly  these  diffracted 
radial  rays  do  not  all  bear  equal  energy,  as  they  did  in  the  unbounded  medium,  but 
have  intensities  appropriate  to  their  positions  in  the  diffraction  pattern. 

The  piston  problem  illustrates  several  interesting  properties  of  exact  ray  the¬ 
ory.  It  demonstrates,  as  did  the  harmonic  oscillator  waveguide  Ccise,  that  fruitful 
ray  analysis  can  sometimes  proceed  even  when  there  are  large  regions  where  the  ray 
paths  cannot  be  computed.  In  exact  ray  theory,  problem  areas  may  sometimes  be 
“skipped  over,”  whereas  in  classical  ray  theory,  ray  traces  must  usually  proceed  from 
the  source  without  interruption. 

The  piston  problem  involved  a  distributed  source,  an  often  problematical 
situation  for  classical  ray  theory  because  it  is  not  always  clear  how  to  initiate  the 
ray  trace.  In  the  exact  ray  theory  analysis  we  simply  stepped  over  the  troublesome 
nearfield  and  took  advantage  of  knowledge  of  the  farfield  results. 

The  method  of  analysis  of  the  piston  problem  is  also  applicable  to  the  problem 
of  a  plane  wave  normally  incident  on  a  circular  hole  in  a  baffle.  If  the  hole  is  not  so 
small  that  Aru  1,  then  the  acoustic  field  in  the  hole  may  be  taken  to  be  that  of  the 
incident  field  (the  Kirchhoff  approximation)  and  the  analysis  proceeds  cis  it  did  for 
the  piston  problem.  According  to  classical  ray  theory,  the  rays  passing  through  the 
hole  do  not  diverge,  but  form  a  beam  of  parallel  rays.  The  acoustic  field  is  uniform 
across  this  beam  and  vanishes  abruptly  at  the  edge  of  the  beam.  By  contrast,  exact 
ray  analysis  accounts  for  diffraction  by  predicting  diverging  rays  of  variable  intensity. 

H.  BEAM  DISPLACEMENT  AND  REFLECTIONS  FROM  A  BOUNDARY 

Until  recently  most  computer  implementations  of  ray  theory  have  treated  acous¬ 
tic  interaction  with  the  surface  and  bottom  boundaries  in  a  rather  simplistic  fashion, 
letting  the  rays  reflect  in  the  specular  direction  and  imposing  a  plane  wave  reflection 
loss  and  phase  shift,  and  possibly  a  frequency  dependent  rough  scattering  loss.  This 
approach  works  quite  well  when  applied  to  a  smooth  ocean  surface,  because  then 
the  surface  is  an  almost  perfect  pressure  release  boundary;  the  rays  reflect  in  the 
specular  direction  with  a  180°  phase  shift  and  undiminished  amplitude.  (Attempts 
to  treat  scattering  from  a  rough  sea  surface  have  had  mixed  success.)  The  method 
meets  with  indifferent  to  poor  success  when  applied  to  bottom  interaction.  It  is  noto¬ 
riously  inadequate  in  shallow  water,  where  bottom  interacting  energy  can  dominate 
the  propagation  problem. 

The  lateral  displacement  of  a  bounded  beam  reflecting  from  an  interface  is  a 
well  kjiown  phenomenon. [33]  In  1981,  Tindle  and  Bold[3'1]  published  their  discovery 
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that  ray  theory  calculations  could  be  made  very  accurate,  even  in  shallow  water,  by 
imposing  a  beam  displacement  and  phase  shift  on  bottom  reflected  rays,  so  that  the 
reflected  ray  emerged  from  the  bottom  somewhat  displaced  from  its  point  of  incidence 
(see  Fig.  4.10),  and  with  an  additional  phase  shift  imposed. 


•  source 


sediment  layer 
interface 


Figure  4.10  Beam  displacement  of  a  bottom  reflected  ray. 


Implementation  of  beam  displacement  in  a  computer  ray  model  is  straightfor¬ 
ward.  The  complex  plane  wave  reflection  coefficient  V^(^)  is  computed  for  a  ray 
incident  on  the  bottom  at  grazing  angle  9  in  the  usual  fashion.  But  the  point  of 
origin  of  the  reflected  ray  is  displaced  from  the  point  of  incidence  by  a  distance  Ar 
given  by 


Ar  = 


1 

kb  sin  9  69  ’ 


(4-34) 


where  kb  is  the  wave  number  in  the  water  at  the  bottom  interface,  and  4>  is  the  phase 
of  V .  There  is  also  a  beam  displacement  phase  shift  A^i,  given  by 

A$(,  =  fefcAr  cos  9 


in  addition  to  the  bottom  reflection  phase  shift  directly  attributable  to  V.  Caustics 
will  usually  be  created  by  including  beam  displacements  even  if  none  were  present 
before,  and  these  must  be  treated  by  uniform  asymptotics.  Curiously,  the  phase  shift 
accrued  by  a  ray  which  grazes  a  beam  displacement  caustic  is  7r/2  (a  phase  had) 
instead  of  the  familiar  — 7r/2  phase  lag  from  encountering  a  conventional  caustic. 

All  of  these  seemingly  peculiar  manipulations  have  been  shown  to  be  mathemat¬ 
ically  well  justified. [33,35]  Batorsky  and  Felscn[36]  showed  that  beam  displacements 


are  inherent  in  the  transformation  from  normal  mode  to  ray  theory  representations 
when  a  bottom  is  present.  And  there  are  compelling  computational  reasons  for  pre¬ 
ferring  a  ray  calculation  to  a  normal  mode  calculation,  even  though  the  normal  mode 
result  is  slightly  more  accurate:  ray  theory  models  can  easily  treat  broadband  sources 
and  sloping  bottoms.  Even  so,  beam  displacement  seems  to  many  to  be  an  ad  hoc 
correction  with  several  counterintuitive  properties,  and  this  has  discouraged  its  use. 

He2iding  the  list  of  disquieting  properties  is  the  fact  that  a  single  beam-displaced 
eigenray  can  emerge  from  the  bottom  and  properly  account  for  the  bottom  interaction 
effects  of  an  ocean  bottom  consisting  of  several  partially  reflecting  sediment  layer 
interfaces.  Experience  with  classical  ray  theory  suggests  that  rays  should  split  at 
the  interfaces  to  account  for  the  partial  reflection  and  transmission  of  energy,  and 
that  this  should  result  in  a  large,  or  even  infinite  number  of  eigenrays  emerging  from 
the  bottom.  This  notion  of  how  the  rays  should  behave  persists  even  though  this 
model  of  bottom  interaction  produces  poor  predictions.  Then  again,  the  computer 
modeler  can  elect  to  use  a  split-ray  approach,  but  with  beam  displacements  imposed 
on  every  partially  reflected  ray,[37]  although  the  split-ray  approach  complicates  the 
modeler’s  task.  The  reader  has  probably  already  begun  to  suspect  that  this  single 
eigenray /multiple  eigenray  quandary  is  another  manifestation  of  the  confounding  of 
cw  and  broadband  phenomena. 

1.  Reflection  from  a  boundary 

If  a  point  source  is  located  above  a  rigid  perfect  reflector  as  shown  in  Fig.  4.11, 
then  the  acoustic  field  is  given  by 

'PCr,  z;  2o)  =  +  'I'r.fl,  (4-35) 


where 

^sphere  —  ^  » 

0? 

II 

1 

and 

R^  =  r^  +  {z-  ^o)^ 

Rl  =  r^  +  {z  +  zo) 

(the  source  strength  has  been  dropped  in  Elq.  [4-35]).  Classical  ray  theory  yields  the 
exact  solution  to  the  Helmholtz  equation  in  this  case,  and  the  two  classical  eigenrays 
are  shown  in  Fig.  4.11.  The  direct  path  eigenray  connects  the  source  at  (0,  zq)  with 
the  observation  point  at  {r,z).  The  bottom  reflected  eigenray  may  be  replaced  by 
an  equivalent  ray  which  connects  the  source  image  at  (0,  —Zq)  with  the  observation 
point. 

If  we  proceed  to  trace  the  exact  eigenray  by  solving  dr  Ida  =  Im  as 

wc  have  in  all  the  previous  sections,  the  result  is  the  wavy  curve  shown  in  Fig.  4.11. 
(This  exact  eigenray,  incidentally,  was  traced  backwards  from  the  observation  point 
to  the  source.)  Two  classical  eigenrays  are  replaced  by  a  single  exact  eigenray  in  the 
usual  fashion,  but  the  result  does  not  seem  particularly  illuminating. 

A  more  fruitful  approach  is  suggested  by  the  representation  of  the  total  field 
in  Eq.  (4-35)  as  a  superposition  of  a  spherical  wave  and  a  reflected  wave  field.  Each 
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Figure  4.11  Reflection  from  a  rigid  boundary.  The  classical  eigenrays  (solid 

STRAIGHT  LINE  SEGMENTS)  MAY  BE  REPLACED  BY  A  SINGLE  EXACT  EIGENRAY  (WAVY  CURVE). 


of  the  contributing  wave  fields  can  be  represented  separately  by  an  exact  eigenray. 
Thus  the  exact  eigenray  representing  'J'.phere  is  the  usual  straight  line  from  source  to 
observation  point.  The  exact  ray  representing  ^rea  connects  the  source  image  point 
with  the  observation  point;  it  is  the  equivalent  of  the  bottom  reflected  ray  of  classical 
ray  theory.  The  exact  and  classical  eigenrays  coincide,  of  course,  because  classical 
ray  theory  is  exact  in  this  case.  The  representation  of  the  acoustic  field  by  two  terms, 
each  with  its  own  exact  ray  representation,  is  the  key  to  establishing  the  connection 
between  exact  ray  theory  and  beam  displacement. 

The  phenomenon  of  beam  displacement  manifests  itself  only  when  the  reflect¬ 
ing  boundary  is  not  a  perfect  reflector.  We  will  treat  the  simplest  such  case;  that 
of  a  point  source  in  a  homogeneous  fluid  half-space  overlying  another  homogeneous 
half-space.  The  total  field  is  again  written  as  the  sum  of  a  spherical  wave  field  and 
a  reflected  field,  and  ^sphere  is  computed  as  before.  To  compute  the  reflected  field 
the  spherical  wave  impinging  on  the  interface  is  expanded  in  terms  of  plane  waves 
(the  expansion  takes  the  form  of  an  integral  over  the  grazing  angles).  The  reflected 
field  is  then  obtained  by  computing  the  reflected  plane  wave  corresponding  to  each 
incident  plane  wave  by  application  of  a  plane  wave  reflection  coefficient.  The  contri¬ 
butions  of  the  reflected  plane  waves  are  integrated  by  a  saddle  point  approximation. 
The  calculation  is  lengthy;  fortunately,  the  details  are  supplied  by  Brekhovskikh,[33) 
pp.  242-255,  and  we  need  only  state  the  results: 

gffcRi 

^reft  =  Q~~n  ’ 
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Q  =  V{0)-iM{B,R,), 

M  =  Ni6)lkRu  (4-36) 

N  =  \[V"{6)-V\9)i!me], 

tanB  =  («  +  -2o)/r, 

which  are  valid  for  kR  ^  1  -C  kRi  and  for  grazing  angles  0  not  closely  approaching 
the  critical  angle  0„  below  which  total  internal  reflection  occurs.  Primes  denote 
differentiation  with  respect  to  the  grazing  angle  0  (Brekhovskikh  works  with  the 
angle  of  incidence,  which  accounts  for  the  differences  between  the  formulas  in  the 
reference  and  those  appearing  here). 

The  logarithmic  derivatives  of  required  for  ray  tracing  are 

(4-37) 

d^^ldz  dQldz,,.,^  ,,dRjdz 

+ 

where 

dQ/dr  =  [irM -{z  +  Zo)ir -iN'/kRi)]/Rl 
dQldz  =  [i{z  +  zo)M-\-r{V' -iN'lkRi)]IR\, 

and 


dRildr  =  r/Ri, 
dRi/dz  =  iz  +  zo)IRi. 

Equations  (4-36)  and  (4-37)  are  valid  regardless  of  the  complexity  of  the  ocean 
bottom. 

In  the  case  of  the  fluid  half-space  bottom,  the  plane  wave  reflection  coeffi¬ 
cient  V{0)  is  given  by 


^  Pi  sing-^(n^  -cos^g)*/^ 

Pi  sin  0  -f  p(n*  —  cos^  ^)’/*  ’ 

where  p  and  c  are  the  density  and  sound  speed  in  the  upper  half-space,  px  and  Ci  are 
the  density  and  sound  speed  in  the  lower  half-space,  and  n  =  c/ci.  In  this  instance, 
V  is  independent  of  frequency.  The  beam  displacement,  computed  according  to 
Eq.  (4-34),  is  zero  unless  the  grazing  angle  is  less  than  the  critical  angle,  given  by 


cos0„  =  n. 


By  letting  m  =  pi/p  and  defining 
o  =  m  sin  0,  ^  =  m  cos 


y  =  (n^  -  cos^P)*/^, 


the  required  derivatives  of  V  and  N  are  found  to  be 

a  +  7 

q(1  -V)  +  7^X1  +  ^)  +  2(^  +  l')V' 

a  +  7 

^(1  -V)  +  7'"(1  +  V)-  3(a  -  yW  -t-  3(;3  +  y')V" 

0  +  7  ’ 

i(V''"-l/"tant?-V"8ec*<?), 

where 

0/3 

m^7’ 

^2  _  ^2  _ 

m^7  ’ 

-7'(4  +  37"/7). 

2.  Frequency  dependent  ray  representation  of  beam  displacement 

Since  the  calculation  of  relies  on  a  saddle  point  approximation  the  re¬ 
sulting  rays  are  not  quite  exact  rays.  Nor  are  they  frequency  dependent  rays,  since 
the  reflection  coefficient  V,  given  by  Eq.  (4-38),  is  frequency  independent.  Rather 
than  engage  in  a  pointless  search  for  a  suitable  qualifying  adjective  I  shall  simply 
continue  to  call  these  exact  rays. 

Figure  4.12(a)  shows  a  point  source  located  100  m  above  a  plane  interface, 
and  an  observation  point  located  50  m  above  the  interface  and  450  m  downrange 
from  the  source.  The  fluid  density  ratio  is  pxjp  —  1.25  and  the  sound  speed  ratio  is 
n  =  c/ci=  1500/1750.  The  solid  straight  line  segments  are  the  classical  direct  path 
and  bottom  reflected  eigenrays.  The  broken  line  is  the  bottom  interacting  eigenray 
with  beam  displacement  included  (the  beam  displacement  was  computed  according 
to  Eq.  [4-34]). 

Figure  4.12(b)  shows  the  frequency  dependent  eigenray  for  which  was 
traced  backwards  from  the  observation  point  until  it  intersected  the  interface.  The 
classical  eigenrays  are  drawn  in  for  comparison.  The  frequency  dependent  eigenray 
for  the  reflected  field  coincides  almost  exactly  with  that  segment  of  the  beam  dis¬ 
placed  eigenray  which  connects  the  bottom  with  the  observation  point,  even  though 
the  methods  used  to  compute  the  exact  ray  and  the  beam  displaced  ray  are  quite 
different. 

The  homogeneous  half-space  bottom  was  chosen  for  the  sake  of  simplicity,  but 
the  substitution  of  a  more  complicated  bottom  would  merely  change  the  details.  Once 
the  decision  is  made  to  represent  the  total  field  as  a  superposition  of  a  waterborne  field 
and  a  bottom  reflected  field  (and  this  is  always  permissible),  then  the  waterborne  and 
reflected  fields  at  an  observation  point  can  always  be  represented  separately  by  their 
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Figure  4.12  Beam  displacement,  (a)  The  solid  straight  line  segments  are  the  clas¬ 
sical  EIGENRAYS.  THE  BROKEN  LINE  IS  THE  BOTTOM  REFLECTED  EIGENRAY  WITH  BEAM  DIS¬ 
PLACEMENT.  (b)  The  frequency  dependent  eigenray  representing  the  reflected  field, 
SHOWN  WITH  the  CLASSICAL  RAYS  FOR  REFERENCE. 
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respective  frequency  dependent  eigenrays.  Moreover,  there  will  always  be  exactly 
one  eigenray  to  account  for  bottom  interaction  in  this  representation  no  matter  how 
complicated  the  bottom  aind,  in  particular,  no  matter  how  many  sediment  layer 
interfaces  there  are.  And  that  eigenray  will  ordinarily  exhibit  a  displacement. 

On  the  other  h2Uid,  it  is  also  permissible  to  write  the  reflected  field  as  a  super¬ 
position  of  reflected  fields  from  each  layer  interface  in  the  bottom.  In  that  case,  each 
of  these  layers  may  be  represented  by  its  own  frequency  dependent  eigenray,  each  with 
its  own  displacement.  Thus  both  approaches  to  implementing  beam  displacement  in 
a  classical  ray  model  are  vindicated. 

It  is  important  to  note  that  the  displacement  of  the  bottom  reflected  frequency 
dependent  eigenray  appeared  naturally;  it  was  not  necessary  to  introduce  it,  as  is  the 
case  in  conventional  ray  tracing. 

I.  EXACT  RAY  THEORY  AND  MULTIPATHINO 

We  saw  in  the  previous  section  that  frequency  dependent  ray  diagrams  could  be 
constructed  for  each  term  in  the  two- term  representation  of  the  field  for  a  source  above 
a  partially  reflecting  boundary.  In  the  case  examined,  the  field  at  an  observation  point 
was  given  by  the  coherent  addition  of  two  frequency  dependent  eigenrays,  one  being 
a  direct  path  ray  from  the  source  to  the  observation  point,  and  the  other  emerging 
from  the  bottom  (see  Fig.  4.12[b]).  Although  the  field  could  also  be  represented  by 
a  single  frequency  dependent  eigenray,  as  in  Fig.  4.11,  it  proved  more  illuminating 
to  construct  two  frequency  dependent  eigenrays.  By  drawing  both  rays  on  the  same 
ray  diagram  we  reintroduce  a  form  of  multipathing. 

In  general,  if  an  expression  for  the  acoustic  field  takes  the  form  of  a  series,  then 
each  term  in  the  series  may  be  represented  by  its  own  ray  diagram.  Such  term-by- 
term  representations  of  the  field  can  be  illuminating,  eis  we  have  just  seen.  If  the 
diagrams  are  so  that  all  of  the  rays  appear  on  the  same  diagram,  then  the 

rays  will  generally  cross  and  a  semblance  of  multipathing  appears.  But  the  fact  that 
a  ray  from  one  diagram  might  cross  or  graze  a  ray  from  another  diagram  should  not 
cause  alarm  that  the  transport  equation  is  poised  once  again  to  predict  singularities 
at  caustics.  Each  ray  diagram  is  actually  independent  of  all  the  others.  No  single 
ray  diagram  will  exhibit  multipathing,  and  the  transport  equation  must  be  applied 
separately  to  each  diagram. 

When  multipathing  is  introduced  in  this  way  the  exact  eigenrays  have  the  peculiar 
property  that  not  all  of  them  can  be  traced  back  to  the  original  source.  The  example 
of  beam  displacement  examined  in  the  last  section  showed  an  exact  ray  apparently 
emerging  from  the  bottom,  but  it  is  quite  meaningless  to  inquire  how  the  ray  “got 
there.”  The  harmonic  oscillator  waveguide  supplies  another  example.  In  that  case, 
the  acoustic  field  was  expressed  as  a  normal  mode  summation.  We  can  now  see  that 
each  term  in  that  series  would  have  its  own  exact  ray  representation  and,  in  fact,  the 
behavior  of  those  rays  is  easily  discerned.  The  rays  for  each  term  in  the  series  are 
simply  horizontal  lines  projecting  from  the  z  axis.  The  2  axis  acts  as  a  line  source 
for  these  rays  even  though  the  actual  source  was  a  point.  We  find  that  K  =  ^'n  for 


the  nth  term,  where  is  the  eigenvalue  for  that  term  and  mode. 


V.  CONCLUSION 


The  computational  convenience  and  physical  insight  provided  by  ray  theory  mod¬ 
els  have  encouraged  their  use  even  under  conditions  where  the  approximation  inherent 
in  the  theory  becomes  a  serious  liability.  In  underwater  acoustics,  for  example,  ray 
models  are  used  increasingly  in  long  range,  low  frequency  propagation  analysis,  and  in 
shallow  water  propagation  analysis.  But  in  these  environments,  ray  theory  deficien¬ 
cies  in  the  treatment  of  caustics,  bottom  interaction,  and  weakly  ducted  propagation 
detract  from  its  usefulness.  In  an  effort  to  retain  the  advantages  of  ray  theory  while 
correcting  its  lapses,  researchers  have  developed  several  wave  theoretical  corrections 
to  classical  ray  theory. 

The  author  has  been  engaged  for  several  years  in  applying  ray  theory  to  the  Helm¬ 
holtz  equation  which  describes  the  acoustic  field  generated  by  a  time  harmonic  source. 
During  the  doctoral  research  which  lead  to  this  dissertation,  I  became  intrigued  by 
the  possibility  that  the  ray  theory  approximation  might  be  abolished  altogether  in 
solutions  to  the  Helmholtz  equation.  The  idea  was  to  recast  the  Helmholtz  equation 
in  a  ray  theoretical  formulation  in  which  a  geometrical  energy  conservation  law  (the 
transport  equation)  would  be  obeyed  exactly.  The  reformulation  in  fact  turned  out 
to  be  very  straightforward,  but  did  not  lead  to  a  practical  computational  scheme  for 
solving  the  Helmholtz  equation.  What  did  emerge,  however,  v/as  a  practical  com¬ 
putational  method  for  constructing  exact  ray  diagrams  from  known  solutions  to  the 
Helmholtz  equation.  Thus,  whenever  a  solution  to  the  Helmholtz  equation  is  avail¬ 
able,  whether  by  analytical  or  numerical  methods,  the  solution  may  be  represented 
in  the  form  of  an  exact  ray  diagram.  The  transport  of  energy  within  the  flux  tubes 
formed  by  exact  ray  paths  is  rigorously  conservative;  the  ratio  of  areas  law,  familiar 
from  classical  ray  theory,  which  governs  intensity  in  infinitesimal  ray  bundles  ceases 
to  be  an  approximation. 

Yet,  in  a  series  of  numerical  experiments  described  in  this  dissertation,  the  exact 
ray  diagrams  constructed  for  several  known  solutions  to  the  Helmholtz  equation  are 
utterly  unlike  the  corresponding  classical  ray  diagrams.  This  is  so  despite  the  fact 
that  classical  ray  theory  makes  quite  accurate  field  predictions  for  most  of  these 
cases,  provided  allowances  are  made  for  caustics  and  other  known  problematical 
regions.  Moreover,  the  differences  are  not  reconciled  at  high  frequencies,  as  might  be 
expected.  The  surprising  fact  is  that,  by  expunging  the  ray  theory  approximation 
altogether,  the  character  of  the  ray  paths  is  changed  completely.  Actually,  several 
principal  features  of  exact  ray  diagrams  had  been  anticipated  long  before  methods 
were  invented  for  tracing  exact  rays,  but  that  hardly  abated  the  sense  of  surprise 
upon  eventual  confirmation. 

The  lack  of  multipathing  is  the  single  feature  which  most  distinguishes  exact  ray 
diagrams  from  their  classical  counterparts.  The  explanation  for  this  difference  goes 
directly  to  the  nature  of  classical  ray  theory  and  the  ray  theory  approximation.  It  is 
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shown  that  classical  ray  theory  is  based  on  the  phenomenology  of  propagating  fronts 
generated  by  pulsed  sources,  while  exact  ray  theory  is  founded  on  the  phenomenology 
of  time  harmonic  fields.  When  classical  ray  theory  is  applied,  somewhat  out  of  its 
natural  context,  to  the  Helmholtz  equation,  it  takes  the  form  of  a  nonuniform  high 
frequency  asymptotic  approximation.  The  nonuniformity  largely  accounts  for  the 
differences  between  classical  and  exact  ray  theory. 

As  to  practical  considerations,  even  though  exact  ray  theory  does  not  presently 
offer  any  improved  computational  methods,  it  does  provide  a  new  method  for  display¬ 
ing  the  acoustic  field.  If  the  new  ray  diagrams  were  merely  minor  perturbations  of 
the  classical  ray  diagrams  then  they  would  be  of  little  interest.  But  because  excict  ray 
diagrams  are  usually  profoundly  different  from  classical  diagrams,  they  can  provide 
insights  offered  by  no  other  method  of  display.  The  computational  method  for  tracing 
exact  rays  lends  itself  to  implementation  in  existing  numerical  propagation  models. 
Suitably  modified  models  would  be  capable  not  only  of  predicting  the  amplitude  and 
phase  of  the  acoustic  field  but  also  of  tracing  out  its  lines  of  propagation. 
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Slidell,  LA  70458 
Attn:  J.  Paquin 

M.  Bradley 
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Attn:  Lan  Nghiem-Phu 

TRW,  Inc. 
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Attn:  R.  T.  Brown 
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Geophysical  and  Polar  Research  Center 
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1215  W.  Dayton  Street 
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Honolulu,  HI  96822 

Attn:  L.  N.  Frazer 
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Director 
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9  Grove  Street 
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CANADA 

Attn:  D.  Chapman 
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Rosenteil  School  of  Marine  and  Atmospheric  Science 
Division  of  Ocean  Engineering 
University  of  Miami 
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Attn:  H.  DeFarrari 

T.  Yamamoto 

Applied  Physics  Laboratory 
The  Johns  Hopkins  University 
John  Hopkins  Road 
Laurel.  MD  20810 
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Massachusetts  Institute  of  Technology 
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The  Catholic  University  of  America 
Washington,  D.  C.  20064 

Attn;  H.  Uberall 

J.  J.  McCoy 
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The  University  of  Miami 

10  Rickenbacker  Causeway 

Miami,  FL  33149 

Attn:  F.  Tappert 
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Department  of  Electrical  Engineering 
Polytechnic  Institute  of  New  York 
Farmingdale,  NY  11735 

Attn;  L.  B.  Felsen 
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1.  Tolstoy 

Knockvennie,  Castle  Douglas 

S.W.  SCOTLAND 

GREAT  BRITAIN 
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Department  of  Aerospace  Engineering 
and  Engineering  Mechanics 

The  University  of  Texas  at  Austin 
Austin,  TX  78712 
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The  University  of  Texas  at  Austin 
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